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ABSTRACT 


The  Naval  Weapons  Support  Center  is  planning  to  implement  a  bonus  system  to 
improve  the  reliability  of  pyrotechnic  devices.  The  measure  of  effectiveness  that  they 
wish  to  use  to  determine  how  to  award  bonuses  is  the  reliability  of  pyrotechnic  devices. 
The  data  available  to  estimate  this  reliability  is  based  on  the  current  sampling  inspection 
plan  in  which  devices  are  tested  in  different  environments.  The  models  which  include 
both  dependence  and  independence  assumptions  between  the  outcomes  of  these  tests  are 
implemented  and  estimates  of  overall  reliability  along  with  95  %  lower  confidence  bound 
are  obtained.  The  95  %  lower  confidence  bounds  are  found  by  bootstrapping.  Using 
these  estimates,  models  for  making  the  decision  to  award  bonuses  are  discussed  and 
studied  using  Monte  Carlo  simulation  . 
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I.  INTRODUCTION 


A.  BACKGROUND 

Nations  spend  a  lot  of  money  to  establish  a  strong  defense  network.  It  is  essential 
that  nations  buy  reliable  weapons  and  ammunition  from  the  contractors.  To  ensure  re¬ 
liability,  contracts  must  include  a  lot  acceptance  sampling  plan  that  specifies  the  minimal 
acceptable  quality. 

The  Naval  Weapons  Support  Center  purchases  pyrotechnic  devices.  Unless  other¬ 
wise  specified  in  the  contract,  the  supplier  is  responsible  to  see  that  his  devices  meet  all 
inspection  requirements  as  specified.  The  inspection  requirements  are  particular  to 
characteristics  of  each  type  of  device  and  are  specified  on  the  reference  drawings  and 
supplemental  quality  assurance  provisions  of  the  contract.  Any  testing  that  needs  to  be 
done  on  these  devices  is  explained  in  this  contract. 

When  nations  buy  weapons  and  ammunition  from  the  same  contractor,  they  would 
like  the  quality  to  improve  over  time.  As  contracts  are  now  written,  contractors  need 
only  to  satisfy  the  requirements  of  the  sampling  inspection  plan  for  lot  acceptance. 
Under  such  contracts,  contractors  have  no  incentive  to  improve  the  quality  of  items  they 
provide.  For  this  reason,  to  improve  quality,  The  Naval  Weapons  Support  Center  has 
decided  to  implement  a  bonus  system.  The  contractor  will  be  awarded  a  bonus  if  the 
result  of  the  sampling  inspection  exceeds  the  minimum  requirements  for  lot  acceptance. 

The  Naval  Weapons  Support  Center  will  begin  to  implement  a  bonus  system  for 
pyrotechnic  devices  in  FY91.  The  data  available  to  make  the  decision  whether  to  award 
a  bonus  is  based  on  the  current  sampling  inspection  plan.  This  plan  is  a  series  of  de¬ 
structive  tests  in  different  environments.  The  purpose  of  this  thesis  is  to  provide  the 
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Naval  Weapons  Support  Center  with  guidance  for  implementing  a  bonus  system  for 
pyrotechnic  devices.  Once  implemented,  the  pyrotechnic  bonus  system  will  serve  as  a 
prototype  for  bonus  systems  for  other  devices. 

B.  PYROTECHNIC  DEVICE  RELIABILITY 

A  pyrotechnic  device  is  a  chemical  and  grenade  ammunition.  There  are  three  cate¬ 
gories  of  pyrotechnic  devices.  The  three  categories  are:  Aerial  display,  Surface  display 
and  Grenades  (Ref.  1:  pp.  I-2J.  Samples  of  the  Pyrotechnic  devices  are  exposed  to 
various  environments  and  then  activated.  The  criteria  for  successful  activation  depends 
on  the  type  of  device  : 

1.  Aerial  display  (Ground  signals,  flares,  airburst  simulators,  sub  signals,  signal  kits, 

etc.) 

"Successful  activation  means  that  the  item  will,  after  simulating  user  environment, 
successfully . " 

•  launch, 

•  have  proper  separation  /  signal  ignition, 

•  reach  desired  altitude  at  correct  angle, 

•  have  proper  parachute  deployment, 

•  have  proper  display  color, 

•  have  proper  display  time, 

•  have  no  subsequent  interference  with  next  item. 

2.  Surface  (Ground  or  Water)  display  (flare,  hand-held  signals,  smoke  /  illume 
grenades,  simulators,  etc.) 

"Successful  activation  means  that  the  item  will,  after  simulating  user  environment, 
successfully 
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•  have  signal  ignition  with  proper  display, 

•  have  proper  display, 

•  have  proper  display  time, 

3.  Grenades  (fragmentation,  defensive,  white  phosphors,  etc.) 

Successful  activation  means  tha*  the  item  will,  after  simulating  user  environment, 
successfully 

•  function  (high  order  detonation  following  delay), 

•  have  proper  dissipation  of  payload, 

•  have  completely  consumed  payload. 

C.  PYROTECHNIC  DEVICE  RELIABILITY  PROBLEM 

Samples  from  any  large  lot  of  pyrotechnic  devices  submitted  by  a  manufacturer  must 
activate  after  exposure  to  different  environments.  These  are; 

1.  Manufacturer  Environment, 

2.  Temperature  and  Humidity  Environment, 

3.  Vibration  Environment, 

4.  Altitude  Environment. 

All  items  tested  are  subjected  to  the  manufacturer  environment.  However  items  are  only 
subjected  to  one  of  the  three  remaining  environments  :  Temperature  and  Humidity, 
Vibration  or  Altitude. 

The  sampling  plan  consists  of  using  four  distinct  samples  from  a  lot  that  can  be  as¬ 
sumed  (approximately)  statistically  independent.  The  items  tested  in  each  sample  are 
also  assumed  to  be  independent.  According  to  the  sampling  plan; 

•  20  items  are  subjected  to  the  Manufacturer  Test, 


3 


•  20  items  are  subjected  to  both  the  Temperature  and  Humidity  Test  and  Man¬ 
ufacturer  Test, 

•  32  items  are  subjected  to  both  the  Vibration  Test  and  Manufacturer  Test, 

•  20  items  are  subjected  to  both  Altitude  Test  and  Manufacturer  Test. 

A  total  of  92  items  are  tested. 

Acceptance  criteria  for  each  test  are  : 

1.  Manufacturer  Test  :  Of  the  20  items;  if  no  more  than  1  fails  to  activate,  the  lot 
passes. 

2.  Joint  Temperature  and  Humidity  and  Manufacturer  Test  :  Of  the  20  items;  if  no 
more  than  1  fails  to  activate,  the  lot  passes. 

3.  Joint  Vibration  and  Manufacturer  Test  :  Of  the  32  items;  if  no  more  than  2  fails 
to  activate,  the  lot  passes. 

4.  Joint  Altitude  and  Manufacturer  Test  :  Of  the  20  items;  if  no  more  than  1  fails  to 
activate,  the  lot  passes. 

The  number  of  failures  for  these  tests  will  be  summarized  by  the  vector  : 

(  FOM,  FOTH,  FOV,  FOA  )  (1.2) 


where 

•  FOM  represents  the  number  of  failures  after  the  manufacturer  test, 

•  FOTH  represents  the  number  of  failures  after  both  the  temperature  and  humidity 
and  manufacturer  tests, 

•  FOV  represents  the  number  of  failures  after  both  the  vibration  and  manufacturer 
tests, 

•  FOA  represents  the  number  of  failures  after  the  altitude  and  manufacturer  tests. 
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For  example 


(  1,  1,  2,  1  ).  (1.2) 

represents,  1  failure  in  manufacturer  test,  1  failure  in  the  joint  temperature  and  humidity 
and  manufacturer  test,  2  failures  in  the  joint  vibration  and  manufacturer  test,  and  1 
failure  in  the  joint  altitude  and  manufacturer  test.  It  is  also  the  maximum  number  of 
failures  in  each  test  that  still  leads  to  lot  acceptance. 

The  marginal  distributions  of  the  number  of  failures  in  each  of  the  above  tests  are 
modeled  by  binomial  distributions.  There  are  24  possible  realizations  of  the  sampling 
inspection;  ranging  from  the  best  case  with  no  failures  to  activate,  to  the  worst  case  with 
I,  1,  2  and  1  devices  failing  to  activate  in  these  tests  manufacturer,  temperature  and 
humidity,  vibration  and  altitude  respectively.  These  cases  are  tabulated  in  Table  1. 

To  aw'ard  bonuses  we  need  one  measure  of  effectiveness  for  pyrotechnic  devices  that 
can  be  estimated  from  the  available  data.  Ideally,  this  measure  is  the  reliability  of  the 
device.  However,  because  each  of  the  tests  are  destructive,  there  is  no  one  natural  defi¬ 
nition  of  reliability  for  these  devices.  To  be  on  the  conservative  side,  we  define  the  reli¬ 
ability  of  a  device  to  be  the  probability  that  the  device  will  activate  after  exposure  to  all 
of  the  environments.  It  is  not  easy  to  estimate  this  reliability  from  the  sampling  plan 
data.  This  data  is  incomplete  in  the  sense  that  we  have  limited  information  about  the 
joint  probability  of  activation  after  exposure  to  more  than  one  environment.  To  try  to 
compensate  for  this  lack  in  the  data,  we  will  use  models  for  the  joint  distribution  of 
(  FOM,  FOTH,  FOV,  FOA  )  that  specify  particular  types  of  dependence  between  the 
events  that  a  device  activates  after  exposure  to  different  environments. 


5 


6 


Using  these  models  and  based  on  sampling  plan  data,  estimates  of  the  overall  reli¬ 
ability  along  with  lower  confidence  bounds  are  obtained.  These  will  be  used  to  imple¬ 
ment  the  bonus  system  for  pyrotechnic  devices.  We  compute  the  maximum  likelihood 
estimator  (  MLE  )  of  the  reliability  by  maximizing  the  appropriate  likelihood;  lower 
confidence  bounds  (  LCB  )  are  found  by  bootstrapping.  The  M  L Li's  are  computed  under 
both  independence  and  dependence  assumptions.  The  estimation  procedures  assuming 
independence  are  described  in  Chapter  II.  In  Chapter  III  we  incorporate  dependence 
by  fitting  a  Log  Linear  Model  to  our  data.  The  MLE's  from  Chapter  II  and  Chapter 
III  lead  to  inappropriate  results  for  this  for  this  problem;  thus  in  Chapter  IV  we  consider 
alternate  and  very  conservative  estimates  of  reliability.  Using  the  estimates  of  Chapter 
IV,  we  investigate  a  sequential  scheme  for  making  the  decision  to  award  bonuses  in 
Chapter  V.  The  results  of  simulations  are  presented  in  Chapter  VI.  Finally,  conclusions 
and  recommendations  are  given  in  Chapter  VII. 
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II.  THE  ESTIMATION  OF  THE  MAXIMUM  LIKELIHOOD  ESTIMATOR 


(MLE)  OF  THE  RELIABILITY  WITH  INDEPENDENCE  ASSUMPTION 
A.  DEFINITIONS 

We  will  say  that  a  device  survives  environment  E  ,  if  it  is  still  potentially  capable  of 
activation  after  exposure  to  environment  E  .  Let, 

•  £,  be  the  device  activates  after  exposure  to  Manufacturer  environment, 

•  £j  be  the  device  activates  after  exposure  to  Temperature  and  Humidity  environ¬ 
ment, 

•  £3  be  the  device  activates  after  exposure  to  Vibration  environment, 

•  £„  be  the  device  activates  after  exposure  to  Altitude  environment. 

We  define  the  reliability  of  device  as  below, 

r  =  p ( £,  n  e,  n  £3  n  £*).  (2.0 

In  this  formula,  R  means  the  probability  that  a  device  activates  after  exposure  to  four 
environments.  We  will  estimate  R  for  each  of  the  24  cases  which  lead  to  lot  acceptance. 

Let  <2.  —  P  (  £, )  be  the  probability  that  device  activates  after  exposure  to  environ¬ 
ment  i,  for  i  =  1,  2,  3,  4.  In  the  acceptance  sampling  plan  several  of  the  items  must 

activate  after  exposure  to  a  joint  manufacturer  and  another  environment.  To  avoid 
confusion  we  will  denote  tests  1  through  4  as  the  manufacturer  test,  the  joint  temper¬ 
ature  humidity  and  manufacturer  test,  the  joint  vibration  and  manufacturer  test  and 
joint  altitude  and  manufacturer  test  respectively. 
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Let 


£,  =  £(£,)  =  {?,  (2.2) 

and  let 

Ri  =  P(E,n  Et)  .  (2.3) 

Here  R ,  is  the  probability  that  a  device  survives  test  i  for  i  =  1,  2,  3,  4.  The  simplest 
model  is  to  assume  that  E £3,  £4  are  independent.  If  we  assume  that  £,....,  £4 
are  independent  then 

Ri  =  <?,  a  (2-4) 

for  i  =  2,  3,  4  and  the  reliability  of  device  is, 

R=  Qx  QiQiQ*-  (2-5) 

B.  THE  LIKELIHOOD  EQUATION 

Let 

•  T,  be  the  number  of  devices  that  activate  after  test  i, 

•  n,  be  the  number  of  items  given  test  i. 

Then  X,  is  binomial  with  parameters  R,  and  n,  for  i=  1,  2,  3,  4.  Under  the  assumption 
of  independence  the  joint  likelihood  function  of  observing  Xl  =  jc,  ,  .  .  .  ,  A4  =  x,  is 


L  {  x]  ,  x2 ,  x3  ,  \  R}  ,  R2 ,  R3 ,  R4) 


i=i 


(2.6) 
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with  constraints  : 


0  <  R]  <  1, 


0  <  R2  <  Ru 


0  <  /?3  <  /?,, 

0  <  /?d  <  R}. 

Our  aim  is  to  maximize  this  likelihood  function  subject  to  the  constraints  that 
(  ,  R2 ,  Rz ,  /^ )  e  S  where 

S  =  { (  R]  R2  R3  ) :  0  <  /?,  <  1,  0  <  /?,  <  /?,  i  =  2,  3,  4  }  . 

From  the  equation  (2.6),  we  see  that  maximizing  L  is  equivalent  to  maximizing 


/  =  £{(*,-  In/?,)  +  (nt  -  Xi)  ln(  1  —  /?,)} ,  (2.7) 

i=l 


,,rhere  the  constant  multipliers 


for  i  =  1,  2,  3,  4  have  been  dropped  (  because  they 


do  not  effect  the  maximization  procedure  )  and  the  natural  logarithm  of  L  is  taken. 

We  first  show  that  /  is  a  concave  function.  To  show  that  /  is  a  concave  function, 
we  can  show  -  /  is  a  convex  function.  According  to  Theorem  3.3.6  [Ref.  2:  p.  92], 
by  looking  at  its  Hessian  matrix,  we  can  learn  whether  function  is  convex  or  not.  If  its 
Hessian  matrix  is  positive  semi-definite  at  each  point  S  then  function  /  is  convex.  To 
create  the  Hessian  matrix,  we  must  calculate  partial  derivatives  of  the  function  -  / , 


vl  __fj_  n<  ~  Xi 

dRt  R,  1  -  Rt 


(2.8) 
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(2.9) 


Al_  =  jl  .  n‘  ~  x‘ 

dR-  Rf  1  -  R* 


Then,  the  determinant  of  the  Hessian  is: 


(2.10) 


4 


Clearly  we  can  see  that  for  0  <  R,  <  1  i  =  1,  2,  3,  4,  |  H  |  is  always  positive.  Because 
-  /  is  continuous,  this  implies  that  -  /  is  a  convex  function  on  S.  As  a  result  of  this, 

/  is  a  concave  function. 

We  note  that  with  the  constraints  on  the  probabilities  Rx  ,  .  . .  ,  Kx  ,  there  does  not 
in  general  exist  a  closed  form  solution  to  MLE.  However,  with  only  24  realizations  of 
xx  ,  .  .  .  ,  xa  of  interest,  the  estimated  reliabilities  for  these  24  cases  can  be  found  with 
some  rather  tedious  but  straight-forward  computations. 

C.  COMPUTING  THE  MAXIMUM  LIKELIHOOD  ESTIMATORS  AND  LOWER 
CONFIDENCE  BOUNDS 

Because  /  is  concave  over  the  convex  set  S,  if  the  maximum  occurs  in  the  interior 
of  S,  it  is  a  unique  maximum  and  is  given  by 

=  4-  (2-12) 

for  i  =  1,2,  3,  4.  where 
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0  <  /?,  <  1  ,  0  <  R'  <  R} 

for  i  =  2,  3,  4. 

In  this  case  the  MLE's  for  Q,  i  =  1,  2,  3,  4  are 

2i  -  %  (2-13) 

iL 

a  n. 

Qi  =  -jcf-  (214) 

~ 

for  i  =  2,  3,  4.  Finally  we  can  estimate  the  reliability  of  the  device  as 

A  A  A  A  A 

R-QiQiQiQ*-  (2-15) 

Table  2  summarizes  the  cases  for  which  the  MLE's  can  be  found  using  (2.12)  - 
(2.15).  When  the  maximum  falls  on  the  boundary  of  S,  there  is  no  explicit  expression 
for  the  MLE  of  (  /?, ,  . . . ,  R* ) .  These  are  the  cases  with  the  exception  of  (  1  1  2  1  ) 
which  have  one  failed  item  in  the  manufacturer  test.  This  implies  that  R ,  will  be  less 
than  1.0  .  To  find  the  MLE,  we  find  the  maximum  of  ( /?,  ,  . . ,  /?,)  on  each  of  the 
boundaries,  compute  /  for  each  of  these  and  let  the  MLE  be  the  one  with  the  largest 
value  of  /  .  It  is  clear  that  in  most  cases  several  of  the  boundaries  can  be  eliminated 
from  consideration,  simplifying  computation  considerably. 
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Table  2.  EASY  AND  HARD  CASES 


EASY  CASES 

HARD  CASES 

( 

0 

0 

0 

0  ) 

(  1 

0 

D 

0  ) 

( 

0 

0 

D 

1  ) 

(  1 

0 

1  ) 

( 

0 

0 

1 

0  ) 

(  1 

0 

1 

0  ) 

( 

D 

1 

0 

0  ) 

(  1 

1 

D 

0  ) 

( 

0 

1 

1  ) 

(  1 

D 

1 

1  ) 

( 

0 

1 

n 

1  ) 

(  1 

1 

D 

1  ) 

( 

D 

1 

i 

0  ) 

(  I 

1 

1 

0  ) 

( 

D 

2 

0  ) 

(  1 

0 

2 

0  ) 

( 

1 

1 

1  ) 

(  1 

1 

1 

1  ) 

( 

0 

0 

2 

1  ) 

(  1 

H 

2 

1  ) 

( 

D 

1 

2 

0  ) 

(  1 

1 

2 

0  ) 

( 

1 

2 

1  ) 

( 

1 

1 

2 

1  ) 

D.  EXAMPLES 
1.  Easy  Case 

In  this  example  our  failure  vector  is. 


(  0  1  2  1  ). 
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According  to  this  failure  vector,  we  can  write  our  likelihood  equation  by  using  equation 
(2.6)  to  get 


L=  R™(  1  -  )°/?’9(  1  -  R2  )'/?330  (  1  -  R2  )2Rl\  1  -  R,  )'  (2.16) 


0  <  Rl  <  1  ,  0  <  R;  <  /?,  , 


for  i  =  2,  3,  4.  From  the  likelihood  above; 


A  20 
=  lk  =  ^  =  1  0000 

a  x,  19 

R>  =  =  25-  =  °-9500 

A  xt,  30 

^ =  -  if  =  0  9375 

A  X4  lO 

^ =  ^  =  15- =  a950°- 


(2.17) 


a  a 

Results  imply  that  all  R's  for  i  =  2,  3,  4  are  between  0.0  and  /?,,  which  means  that 
constraints  are  met.  Now  we  can  estimate  Q,'s  using  equations  (2.13)  and  (2.14). 


Qi  = 

1.0000 

Qi  - 

0.9500 

1 

=  0.9500 

Qz  - 

0.9375 

1 

=  0.9375 

(2.18) 

&  = 

0.9500 

1 

=  0.9500. 

And  finally,  reliability  of  the  device  can  be  estimated  by  using  equation  (2.15) 
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R  =  (  1.0000)  (0.9500)  (0.9375  )  (0.9500)  =  0.84609375. 


(2.19) 


2.  Hard  Case 

In  this  example  our  failure  vector  is 


(  1  0  2  1  ) 


According  to  this  failure  vector,  we  can  write  our  likelihood  equation  by  using  equation 
(2.6)  as 


L  =  tf,19(  1  -/?,)'  R22°(  1  -  R2  )°  Rf  (  1  -  R2  )2  Rl9(  1  -  /?4  )'  (2.20) 


0  <  R]  <  1  ,  0  <  R,  <  /?, 


for  i  =  2,  3,  4.  From  the  likelihood  above; 


*1  19 

—  =  —  =  0  9500 

/i,  20  v.yjw 

xi  20 

=  10000 


f3 

n3 


30 

32 


0.9375 


ft  =  11 

«4  20 


0.9500. 


(2.21) 


X  X 

As  vou  see  from  the  above  >  — --  thus,  the  MLE  does  not  lie  in  the  interior  of  S. 

n2  n\ 

We  begin  by  considering  the  boundary 
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R  j  —  Ri  —  Rj  —  R4  ■ 


The  likelihood  equation  on  this  boundary  is, 


L  =  Rf(  1  -  Rx  )J 
0  <  /?,  <  1  . 


From  the  likelihood  equation, 


A  A  A  A  fifi 

/?,  =  =  R3  =  R*  =  =  0.9565 . 


Then  value  of  the  likelihood  with  these  estimated  /?,'s  is 


L  =  ( 


88 

92 


-.88 


1  92  } 


0.7*1 


10 


-7 


For  the  boundary 


R\  —  R2  ~  R3  ’  R4  —  R\  < 


the  likelihood  equation  is, 


L  =  <(  1  -  /?,  )3*i9(  1  -  *4  )’ 
0  </?,<!,  0  <  /?4  <  /?,  . 


From  the  likelihood  equation, 


A  A  A  fyQ  A  1  Q 

/?,  =  R2=  r3  =  =  0.9583  /?4  =  ~  «  0.9500 , 


(2.22) 


(2.23) 


(2.24) 
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and  the  value  of  the  likelihood  with  these  estimated  R,'s  is 


=  (  -fr  )69(  JT  )3(  )’9(  )'  =  0.723  x  10"7 


For  the  boundarv 


R,  =  R2  =  R,  ,  R3  <  Rx 


the  likelihood  equation  is, 

L  =  Rx*  (  1  -  Rx  f  R™  (  1  -  R3  f 
0  <  R]  <  1  0  <  R3  <RX. 

From  the  likelihood  equation, 

AAA  ZQ 

/?,  =  R2  =  R4  =  -g-  =  0.9666 
A  10 

-  0.9375 . 

Then  value  of  the  likelihood  with  these  estimated  R's  is 

1  =  ( io )58(  ~k )2(  12 )30(  4" )2  =  °-876  x  10-7 ’ 

For  the  boundary 

R\  —  R3  =  R<x  ,  R2  —  < 


the  likelihood  equation  is, 


(2.25) 


(2.26) 


(2.27) 
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(2.28) 


L  =  R\&  (  1  -  Rx  f  R  f 

0  <  Rx  <  1  0  <  R2  <  Ri  . 

Clearly  L  is  maximized  on  the  boundary  of  the  constraints  (2.28),  i.e.  on  another  of  the 
boundaries  of  S  thus  we  can  eleminate  this  case  from  consideration.  For  the  boundary 

R\  —  /?2  i  R3  —  R\  i  R4  —  R\ ) 


the  likelihood  equation  is, 


L  =  Rx9  (  1  —  Rx  )'  R?  {  \  -  R3f  R\\  \  -  R4  )'  (2.29) 


0  <  Rt  <  1  ,  0  <  R3  <RX  ,  0  <  R4  <  R}  . 


From  the  likelihood  equation, 

A  A  -IQ  A  -If)  AlQ 

R  =  R  ~  Ji.  =  0.9750  R3  =  —  =  0.9375  R^  =  =  0.9500 . 

1  40  J  32  20 


Then  value  of  the  likelihood  with  these  estimated  R's  is 


L  -  (  -35-  )”(  )'(  -§■  )!°(  J2  )!(  )'’(  -|r  >'  -  °"°  x  I0‘7<2.30) 


Quick  inspection  reveals  that  the  remaining  boundaries  can  be  eliminated  from  consid¬ 
eration. 

After  boundary  analysis,  we  can  see  that  equation  (2.30)  gives  us  the  maximum 
likelihood  value.  Thus  the  MLE’s  for  /?, ,  R2 ,  Rz ,  R^  are 

Rx  =  0.9750 


18 


R2  =  0.9750 
R2  =  0.9375 

=  0.9500 . 


Now  we  can  estimate  Q' s  with  equations  (2.13)  and  (2.14) 


£>i  =  0.9750 


ft 

ft 

ft 


0.9750 

0.9750 

0.9375 

0.9750 

0.9500 

0.9750 


1.0000 
0.9615 
0.9744  . 


(2.31) 


And  finally  reliability  of  the  device  can  be  estimated  by  using  equation  (2.15) 


R  =  ( 0.9750  )(  1.0000)  (0.9615)  (0.9744)  =  0.91346100.  (2.32) 


We  computed  MLE's  with  a  FORTRAN  program  MLEA  given  in  Appendix 
A.  After  24  replications  of  this  program,  we  can  estimate  R,  s  for  i  =  1,  2,  3,  4  for  each 
of  the  24  outcomes  of  the  sampling  plan  that  lead  to  lot  acceptance.  These  are  given  in 
Table  3,  R  is  given  in  Table  4. 

We  compute  lower  confidence  bounds  for  each  possible  case  by  bootstrapping. 
The  bootstrap  can  be  used  to  produce  approximate  confidence  intervals  in  an  automatic 
way.  There  are  several  ways  to  set  approximate  confidence  intervals  with  bootstrapping. 
These  are  the  percentile  method,  the  standard  method,  the  bias-corrected  percentile 
method  and  the  nonparametric  method  [Ref.  3  :  pp.  67-70).  Let  6  be  an  unknown  pa- 

A  A 

rameter  with  estimator  9.  To  bootstrap,  samples  are  generated  using  6  in  place  of  the 
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unknown  parameter  0.  For  each  of  these  samples  an  estimator  O'  is  computed.  We 
define  G  (  s )  to  be  the  parametric  bootstrap  cumulative  distribution  function  of  O',  i.e. 

A  A 

G  ( s )  is  the  empirical  distribution  of  the  0"s.  All  methods  mentioned  above  use 
percentiles  of  G  to  define  the  confidence  interval.  They  differ  in  which  percentiles  are 
used.  The  percentile  method  was  used  in  our  calculations. 

If  we  use  the  notation  0(a)  for  the  level  a  endpoint  of  one  sided  lower  confi¬ 
dence  interval  for  6  ,  then 

P[0(a)<0]  =  ot.  (2.33) 

An  estimate  of  6  (  a  )  from  the  bootstrap  cumulative  distribution  is  given  by 

0(a)  =  G_i  (a).  (2.34) 

To  get  LCB's  for  R,  for  each  of  the  24  realizations  of  the  sampling  plan  ,  we 
generate  bootstrap  samples  of  random  failure  vectors,  in  which  failures  come  from  in¬ 
dependent  binomial  distribution  with  parameters  (  n, ,  R, )  for  i  =  1,  2,  3,  4.  This  is  done 
using  the  FORTRAN  program  RANVEC  (Ref.  4]  in  Appendix  B.  We  generate  5000 
failure  vectors  for  each  case.  And  then  we  estimate  R  s  from  each  of  the  5000  failure 
vectors  by  the  means  of  the  program  MLEA.  The  next  step  is  to  compute  the  order 

A  A  A 

statistics  of  R 's  from  i?,  to  R5oao.  We  get  the  parametric  bootstrap  cumulative  distrib- 

A  A 

ution  function  (  i.e.  the  empricai  distribution  of  Rx ,  . . . ,  )  with  this  computation. 

This  is  done  using  the  FORTRAN  program  SORT  in  Appendix  C. 

Finally  we  can  get  the  95  %  lower  confidence  bound  using  equation  (2.34)  from 
this  routine.  Reliabilities  and  95  %  lower  confidence  bounds  are  listed  in  Table  4.  Re- 

A 

suits  in  Table  4  are  given  in  descending  order  and  R's  and  LCB's  are  not  ordered  as  we 
expected  them  to  be.  For  example,  the  failure  vector  that  has  the  maximum  number  of 


failure  in  each  test  is  in  the  middle  of  the  table  with  respect  to  R.  It  has  also  a  bigger 
95  %  LCB  than  the  failure  vector  that  has  a  total  of  3  failed  items  in  each  of  the  joint 
tests.  These  results  are  counter-intuitive  because  we  expect  that  more  failures  indicate 
a  lower  overall  reliability.  This  is  reasonable  because  it  is  likely  that  a  device  that  is 
poorly  constructed  is  more  likely  to  fail  in  any  of  the  four  environments.  From  the  re¬ 
sults  in  Table  4,  it  is  clear  that  an  attempt  to  model  dependence  must  be  made  in  order 
to  get  believable  estimates  of  R. 
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Table  3.  TEST 


CASE 


(0  0  0  0) 


(0  0  01) 


(0010) 


PROBABILITIES  (  A  ) 


1. 0000000 


1.0000000 


1.0000000 


1.0000000 


1.0000000 


1.0000000 


1.0000000 


1.0000000 


0.9687500 


1.0000000 


0.9500000 


(0100) 

1.0000000 

0.9500000 

1.0000000 

1.0000000 

(  10  0  0  ) 

0.9891304 

!■ 

.9891304 

0.9891304 

0.9891304 

(0011) 

1.0000000 

1.0000000 

0.9687500 

0.9500000 

(0101) 

1.0000000 

0.9500000 

1.0000000 

0.9500000 

(01 10) 

1.0000000 

■ 

.9500000 

1.9687500 

1.0000000 

(10  0  1) 

■ 

1.9861111 

■ 

.9861111 

0.9861111 

0.9500000 

(1010) 

0.9833333 

■ 

.9833333 

■ 

.9687500 

■ 

.9833333 

(1100) 

0.9861111 

0.9500000 

■ 

.9861111 

0 

.9861111 

(0020) 

1.0000000 

1.0000000 

■ 

.9375000 

1.0000000 

(0111) 

1.0000000 

■ 

.9500000 

■ 

.9687500 

0.9500000 

(1011) 

0.9750000 

■ 

.9750000 

■ 

.9687500 

■ 

.9500000 

(110  1) 

■ 

.9807692 

0.9500000 

0.9807692 

0.9500000 

(1110) 

0.9750000 

■ 

.9500000 

■ 

.9687500 

.9750000 

(  0  0  21) 

1.0000000 

1.0000000 

■ 

.9375000 

.9500000 

(0120) 

1.0000000 

■ 

.9500000 

■ 

.9375000 

1.0000000 

( 1020) 

■ 

.9833333 

■ 

.9833333 

■ 

.9375000 

.9833333 

(1111) 

■ 

.9615384 

■ 

.9500000 

■ 

.9615384 

.9500000 

(0121) 

1.0000000 

■ 

.9500000 

0.9375000 

.9500000 

(10  2  1) 

0.9750000 

0.9750000 

■ 

.9375000 

.9500000 

(112  0) 

■ 

.9750000 

.9500000 

■ 

.9375000 

.9750000 

(112  1) 

■ 

.9500000 

.9500000 

■ 

.9375000 

.9500000 
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Table  4.  RELIABILITIES  AND  95  %  LOWER  CONFIDENCE  BOUNDS 


(A) 


FAILURE 

VECTOR 

A 

R 

FAILURE 

VECTOR 

95  %  LCB 

(  0  0  0  0  ) 

1.0000000 

(  0  0  0  0  ) 

1.0000000 

( 1000) 

0.9891304 

(0010) 

0.9062500 

(0010) 

0.9687500 

(1000) 

0.9000000 

(1010) 

0.9687499 

1 

1 

0.8750000 

(  0  0  01) 

0.9500000 

(1010) 

0.8550000 

(0100) 

0.9500000 

(  0  0  01) 

0.8500000 

(10  0  1) 

0.9499999 

(0100) 

0.8500000 

(1 100) 

0.9499999 

(10  0  1) 

0.8282812 

(10  11) 

0.9439102 

(1  10  0) 

0.8258822 

(1110) 

0.9439102 

(01 10) 

0.8234775 

(1111) 

0.9385999 

(0011) 

0.8234375 

(  0  0  2  0  ) 

0.9375000 

(1020) 

0.8125000 

(  1  0  2  0  ) 

0.9374999 

(1011) 

0.8015624 

(112  1) 

0.937499 

(0101) 

0.8000000 

■BISinii 

0.9203125 

(1110) 

0.8000000 

(01 10) 

0.9203125 

(0120) 

0.7875000 

(1101) 

0.9201959 

(  0  0  2  1  ) 

0.7749999 

(10  2  1) 

0.9134614 

(1101) 

0.7749018 

(1120) 

0.9134614 

(1111) 

0.7649999 

(0101) 

0.9025000 

(1021) 

0.7647058 

(0021) 

0.8906249 

(1120) 

0.7614843 

(0120) 

0.8906249 

(1121) 

0.7505192 

(0111) 

0.8742968 

(0111) 

0.7499999 

(0121) 

0.8460937 

(0121) 

0.7171874 
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III.  LOG  LINEAR  MODEL  WITH  DEPENDENCE  ASSUMPTION 


A.  BACKGROUND 

Modeling  the  outcomes  of  the  various  environmental  tests  as  independent  is  clearly 
inappropriate.  However,  the  nature  of  the  acceptance  sampling  plan  makes  it  impossible 
to  estimate  the  reliability  of  the  device  (  i.e.  the  probability  that  it  would  activate  after 
exposure  to  all  four  environments  )  without  some  assumptions  about  the  dependence 
between  outcomes  of  various  tests.  Thus,  our  second  approach  is  to  model  the 
pyrotechnic  device  reliability  with  a  log  linear  model. 

The  results  of  our  test  series  create  a(2x2x2x2)  contingency  table.  Let  p 
represent  passing  a  test  and  f  represent  failing  a  test. 


Table  5.  CONTINGENCY  TABLE  STRUCTURE  FOR  TEST  SERIES 


Passed 

Manufacturer  Test 

Failed 

Manufacturer  Test 

Passed 

Tem.&Hum. 

Test 

Failed 

Tem.&Hum. 

Test 

Passed 

Tem.&Hum. 

Test 

Failed 

Tem.&Hum. 

Test 

Passed 

Vibration 

Test 

Passed 

Altitude 

Test 

Mpppp 

MPfPP 

Mfppp 

Mffpp 

Failed 

Altitude 

Test 

MPrrf 

Mpfpf 

Mfppf 

Mffpf 

Failed 

Vibration 

Test 

Passed 

Altitude 

Test 

Mpp/p 

Mpffp 

Mfpfp 

MJJTp 

Failed 

Altitude 

Test 

MppST 

MpJJT 

Mfpff 
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The  results  of  our  tests  can  be  thought  of  as  censored  data  from  a  hypothetical 
(2X2  X  2X2)  contingency  table  (See  Table  1) .  The  frequency  in  each  cell  of  this  table 
is 

U,Jk!,  d,j,  k,  I)  e  {p,f}\ 

the  number  of  devices  out  of  92  which  would  have  result  i  in  environment  1  (manufac¬ 
turer),  result  j  in  environment  2  (temperature  and  humidity  alone),  result  k  in  environ¬ 
ment  3  (vibration  alone),  result  /  in  environment  4  (altitude  alone). 

This  is  a  hypothetical  table,  because  if  a  device  was  exposed  to  all  four  environments 
and  then  failed  to  activate,  there  would  be  no  way  to  discern  which  combination  of  the 
four  environments  caused  failure.  The  data  from  the  acceptance  sampling  plan  can  be 
thought  of  as  censored  data  from  such  a  (2X2X2X2)  contingency  table.  As  an 
example,  a  device  that  is  given  just  the  manufacturer  test  belongs  in  one  of  the  cells 
(  p,  j,  k,  / )  where  (  j,  k,  / )  e  {p,  f}\  because  it  is  not  clear  what  would  have  happened 
to  it  had  it  been  exposed  to  the  other  three  environments. 

Using  a  log  linear  model  under  independence,  the  expected  value  of  each  cell  fre¬ 
quency  is 

ElMyul  =  eM  +  +  ** +  ** +  ^  (3.1) 

where  i.  j,  k,  /  =  p,  f  and 

4  =  -  4  (3.2) 

for  i  =  1,2,  3,  4.  Here  the  parameter  p  represents  the  overall  effect  and  the  parameters 
):e  represent  the  effect  of  passing  environment  i.  It  is  simple  to  show  that  the  log  linear 
model  (3.1)  is  equivalent  to  the  independence  of  the  environments  [Ref.  5:  pp.  25-46). 
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On  a  logarithmic  scale,  the  independence  relation  is  equivalent  to  the  additive  re¬ 


lationship.  As  an  example, 

l°g  Mpppp  =  4  +  '‘■p  +  ^2p  +  J-p  +  ^-p  ■  (3-3) 


Log  linear  models,  which  take  into  account  dependence  include  extra  interaction  terms. 
Because  of  the  extreme  amount  of  censoring,  we  will  only  consider  models  with  two  way 
interaction  terms  of  the  type. 


El  A1Ukll 


.  J1  J2  l3  ,  j4  »  1>2  ,  ,13  ,  ,14  .23  ,  .24  ,  .34 

+  «tj  +  Kk  +  *i/  +  'v*  +  ^  +  ^A/ 


(3.4) 


where  i,  j,  k,  /  =  p,  f, 


(3.5) 


for  i  =  1,  2,  3,  h  and 

4-4= -4- -4 

for  i,  j,  =  1,2,  3,  4.  As  an  example; 


E  [  Mpppp  ]  =  eM  +  kp  +  Xp  +  Xf  +  Xp  +  +  +  ^  +  k”  +  /i"  +  <l" 


(3.6) 


(3-7) 


The  censored  data  does  not  tell  us  much  about  partial  association  between  pairs  of 
tests.  For  this  reason,  we  assume  that  partial  associations  between  any  pair  of  tests  are 
all  the  same.  Let  6  be  the  partial  association  between  tests.  Then  we  can  reformulate 
equation  (3.7)  as  follows: 


El  Mpppp  ]  =  ei‘  +  i\  +  Ji  + 


xl  +  A*  +  (>e 

P  P 


(3.8) 
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The  MLE's  for  the  expected  number  of  devices  in  each  cell  (  or  equivalently  the 
MLE's  for  the  parameters  )  can  not  be  found  explicitly.  The  EM  Algorithm  (  Expecta¬ 
tion  -  Maximization  )  will  be  used  to  approximate  the  MLE's  of  the  expected  number 
of  devices  in  each  cell. 

B.  EXPECTATION  MAXIMIZATION  (  EM  )  ALGORITHM 

A  general  method  of  maximum  likelihood  estimation  from  incomplete  data  is  the 
EM  Algorithm.  The  expectation  maximization  algorithm  is  an  iterative  procedure  where 
each  stage  consists  of: 

■  an  expectation  step  (  E  )  followed  by 
*  a  maximization  step  (  M  ). 

This  algorithm  is  generally  used  to  compute  maximum  likelihood  estimators  in  in¬ 
complete  data  problems.  In  the  application  of  the  EM  algorithm  we: 

•  replace  missing  values  by  their  estimated  expected  values  given  the  incomplete  data 

•  estimate  parameters 

•  reestimate  the  missing  values  assuming  the  new  parameter  estimates  are  correct 

•  reestimate  parameters 

and  so  forth,  iterating  until  convergence  [Ref.  6:  pp.  127-141).  This  iterative  algorithm 
works  as  follows  in  pyrotechnic  device  problem. 

I.  Initialization 

The  algorithm  requires  initial  guesses  for  the  parameters  of  the  log  linear  model 
cell  frequency  in  contingency  table.  It  uses  initial  guesses  and  calculates  initial  cell 
probabilities.  We  will  get  our  initial  guesses  by  first  fitting  a  log  linear  model  under  in¬ 
dependence. 
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2.  Iterations 


=  Expectation  Step 

EM  Algorithm  estimates  expected  cell  frequencies  given  the  data  by  using  the 
most  current  estimates  of  the  parameters.  It  compares  these  estimated  expected  cell 
frequencies  with  the  previous  estimates.  If  the  differences  between  these  two  estimates 
are  small  enough,  then  the  EM  Algorithm  has  converged.  We  then  accept  the  final  es¬ 
timates  of  cell  frequencies  as  MLE,  so  we  can  easily  calculate  cell  probabilities. 

■  Maximization  Step 

If  the  algorithm  has  not  converged,  then  it  starts  to  estimate  new*  parameters  for 
the  log  linear  model  using  expected  cell  frequencies  from  the  E  step  as  if  there  were  ac¬ 
tual  data  available.  Estimation  is  done  by  maximizing  the  likelihood  using  an  iterative 
Newton- Raphson  method.  The  estimated  parameters  from  this  step  are  than  used  in  the 
next  E  step  of  the  EM  algorithm. 

C.  CALCULATIONS 

We  apply  the  EM  algorithm  for  each  realization  of  the  failure  vector  (  FOM, 
FOTH,  FOV,  FOA  ).  During  the  maximization  step  of  the  EM  Algorithm,  we  use  the 
Newton-Raphson  procedure  which  is  described  by  SAS,  [Ref.  7:  pp.  190-212].  We  first 
describe  the  M  step  using  as  an  example  the  following  failure  vector. 

(1  1  2  0) 

■  Maximization  Step 

We  start  with  initial  guesses  (IG)  for  These  are  the  conditional  expected 

values  provided  by  the  previous  E  step.  For  example,  initial  guesses  for  the  failure  vector 
above,  are  shown  below. 
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1G 


M, 


pppp 


XI, 


pppf 


M 


ppfp 


M 


ppff 


M, 


prpp 


M, 


pfpf 


M 


pjjp 


M, 


Pfff 


XI, 


Jppp 


\t , 


fppf 


V, 


fpfp 


M, 


fpff 


M 


ffpp 


Xf, 


'ffpf 


XI, 


Vf/P 


Xf, 


7/// 


83.08 

0.1 

3.36 

0.1 

2.21 

0.1 

0.1 


0.1 

0.1 

0.1 

0.1 

0.1 

0.1 

0.1 


(3.9) 


From  these,  we  can  compute  the  proportion  of  observations  which  fall  into  each  cell  as: 


P 


M, 


pppp 


92 


M, 


pppf 


92 


M, 


ppfp 


92 


%// 

92 


0.9030 
0.00  li 
0.0365 


0.0011 


(3.10) 


We  will  use  P  to  get  MLE's  for  the  parameters  of  the  log  linear  model.  Our  log  linear 
model  which  includes  two  way  interactions  can  be  written  as 
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log  Ppppp 
log  Ppppf 
log  Pppfp 

'll  1  1  16" 

11  1  1-10 

V  ~  P 

'  I  ’ 

1 

1 

log  Pffff 

1  -1  -1  -1  -1  6 

1 

where  yr  =  ( >u,  Xlp,  X2P,  ).],  )*)  is  the  parameter  vector  and  p  is  normalizing  constant 
required  by  the  restriction  that  the  probabilities  sum  to  1  and  where  P,jkl  is  the  proba¬ 
bility  corresponding  to  the  ( / j  k  /)rt  cell  in  the  hypothetical  contingency  table.  We  can 
derive  from  equation  (3.11) 


log  Ppppp 
log  Ppppf 
log  Pppfp 

1  1  1  16' 

1  1  1-10 

/?  -  5 

"  1  " 

1 

1 

Pffff 

-1  -1  -1  -1  6 

1 

where  PT  =  ( X'„  X],  X],  X * )  is  the  parameter  vector  and  <5  is  normalizing  constant  re¬ 
quired  by  the  restriction  that  the  probabilities  sum  to  I.  Then,  to  use  the  SAS  proce¬ 
dure,  we  rewrite  log  Ppppp,  log  Ppppf,  ....  log  Pm  as  the  15  logits,  Fx  =  ( log  PppJ  log  Pmf), 
■  ■  ■  ,  Fx j  =  ( log  PfJ  log  Pm)  so  that 


F(P) 


'  Fl  " 

1  0  0  0  -1 
0100-1 

P’5 

0000-1 

Fpppp 
l°g  Ppppf 

log  Pppfp 

log  Pffff 


(3.13) 
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where  P  =  (  Ppppp,  Pppp/,  . . . ,  Pfm).  Using  equations  above  the  following  result  is  ob¬ 
tained 


F(P) 


2  2  2  2  0 
2220-6 


(3.1-11 


0002-6 


The  design  matrix  (X)  of  this  problem,  from  the  equation  (3.14)  is  as  follows: 


2  2  2  2  0 
2220-6 
2202-6 
2200-8 
2022-6 
2020-8 
2002-8 
2000-6 
0222-6 
0220-8 
0202-8 
0200-6 
0022-8 
0020-6 
0002-6 


(3.15) 


In  the  application  of  Newton-Raphson  Method,  we  use  the  variance  and  covariance 

A 

matrix  S  of  F(  P).  It's  inverse  is  given  by 


S~'(P)  = 


P  -  P 
pppp  1  pppp 


—  P  -X  P 

pppp  1  pppf 

p  _  p2 

*  nnnf  *  nnn/* 


PPPf  rpPPf 


Ppppp  x  Pfffp 

Ppppf  x  Pfffp 


(3.16) 


Ppppp  x  Pfffp 


Ppppf  x  Pfffp 


Pm,  -  P  , 
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The  first  estimates  of  parameters  for  the  Newton- Raphson  procedure  are  calculated 


as  follow's: 


bo  =  lXTS~'{P)XY'  CA'rS“'  (P)F(P)1  (3-17) 

where  S'1  and  F  are  estimated  using  the  proportions  in  the  vector  P.  We  estimate  the 
reduced  logit  response  functions  using  the  equation 

h  =  Xb0.  (3.18) 

From  F0  we  can  compute  the  updated  estimates  of  Ppppp ,  Ppppf ,  ,  Pm  Let 

n  =  (  n,  (  1  ),  n2(  1  ),  . . . ,  ni6  (  1  ) )  be  the  vector,  which  contains  these  estimates  of 

the  cell  probabilities.  The  value  of  the  log  likelihood  evaluated  at  FI  is 


16 

LHE  =  Yj  *i  lQg  n  ( 0 

<=i 


(3.19) 


where  x ,  is  the  number  of  items  in  cell  i  with  respect  to  probability  in  n. 

We  estimate  parameters  iteratively  until  the  difference  between  last  estimate  and 
previous  is  small  enough.  At  each  iteration  we  update  the  inverse  of  variance  and 
covariance  matrix  with  probabilities  of  n  from  the  previous  iteration  and  we  do  fol¬ 
lowing  matrix  computations. 


C  =  J^S-1  (n)* 

G  =  XT\_  92.0  x  (P  -  fl)  ] 


(3.20) 


Let  b,  be  the  next  estimate  of  parameters  /?  in  the  i  th  iteration.  Then  b,  as  follows: 

bt  =  6,_,  -  6  C~'  G  (3.21) 
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where  S  <  1  is  a  constant  supplied  by  the  user. 

We  get  first  cell  expectations  from  initialization.  The  failure  vector,  which  is  defined 
at  page  5,  has  initial  guesses  for  cells  which  are  shown  in  equation  (3.9).  After  the  first 
application  of  Newton  -Raphson  Algorithm,  we  have  following  initial  cell  expectations 
(  CE  ). 


82.4278259 

0.7278117 

3.1047230 

0.1208856 

2.0790358 

0.0809494 

0.3453168 

0.0592893 

2.0268116 

0.0789160 

0.3366428 

0.0578000 

0.2254283 

0.0387050 

0.1651089 

0.1250073 


■  Expectation  Step 

The  conditional  expected  frequency  for  each  cell  expectations  are  calculated  using 
the  estimated  parameters  of  the  log  linear  model  from  the  previous  M  step.  Some  of  the 
cell  expectations  formulated  below 

ElMpm  I  {POM,  FOTH,  FOV,  FOA)  J  «  (20  -  FOM)  x  P  (p  p  p  p  \  p  . . .) 

+  (20  -  FOTH)x  P(pppp\pp  . .) 

+  (32  -  FO V)xP(pppp\p.p.) 

+  (20  —  FOA)  x  P  {p  p  p  p\  p  . .  p) 


ElMpfpf\(FOM,FOTH,FOV,FOA)l  =  (20  -  FO.Xf)  x  P  ( pfpf\  p  .  .  .) 

+  (FOTH)xP(pfpf\pf..) 

+  (32  -  FO  V)  x  P  ( pfpf\  p  ./.) 
+  (FOA)xP(pfpf\p..J) 


where  "  .  "  represents  either  a  pass  or  fail.  For  example, 

P{PPPP\P---)^  the  probability  that  device  passes  all  environments  given  that  it 
passed  manufacturer  test. 

P  (p  P  P  P\  P  P  ■  ^  the  probability  that  a  device  passes  all  environments  given  that 

it  passed  manufacturer  test  and  temperature-humidity  environment. 

The  conditional  probabilities  are  computed  from  the  estimated  cell  probabilities 
from  the  previous  M  step.  As  an  example,  we  have  following  failure  vector, 

(1  1  2  0) 

Then  E  [  Mmr  |  (  1,  1,  2,  0 )  ]  is  calculated  as  follows: 

n,  ,  3  82.4278259 

P(pppp\p...)-  75  -  0.9267193 

,  ,  82.4278259 

KPPPP\PP.-)  -  86.3812462  =  °'9542328 

,  ,  82.4278259  _ _ _ _ 

P(pppp\p.p.)-  ~85~.3I36228  =  09661516 

n/  i  3  82.4278259 

P(pppp\p..p)-  g7-95690|5  -  0.9371388. 


where  the  numbers  come  from  (,->.22)  the  estimated  cell  expected  values  after  from  the 
previous  M  step.  Thus, 
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E  [  Mpppp  1(1,  1,  2,  0)]  =  83.4654163 


The  remaining  16  conditional  expectations  are  calculated  similarly.  Then  we  com¬ 
pare  these  expectations  with  previous  expectations.  If  the  difference  between  compared 
expectations  are  small  enough  then  EM  Algorithm  is  assumed  to  have  converged. 

When  EM  Algorithm  converges,  we  can  estimate  reliability  of  device  (  R  )  after 
exposure  all  four  environments  as 


R  = 


E  [  Mpppp 
92.0 


(3.25) 


For  example,  for  the  failure  vector  above,  after  27  iterations,  we  obtain 

E  [  Mpppp  ]  =  84.412714, 


and 


R  = 


84.412714 

92.0 


0.9175259. 


D.  INITIAL  GUESS  PROBLEM 

The  EM  Algorithm  uses  an  initial  guess  vector  as  in  equation  (3.9).  We  use  the 
independence  assumption  as  described  in  equations  (3.1),  (3.2)  and  (3.3).  The  following 
procedures  are  used  during  the  application  of  EM  Algorithm  to  the  log  linear  model  with 
the  independence  assumption. 

We  have  still  an  initial  guesses  problem  for  parameters  in  the  log  linear  model  with 
the  independence  assumption.  For  this  reason,  we  do  some  calculations  to  find  initial 
guesses  for  the  parameters.  Moreover,  one  can  choose  random  initial  guesses  for  the 
parameters  too. 
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Assume  that  we  have  the  following  failure  vector, 


(10  10). 

In  this  procedure  we  divide  the  number  of  devices  which  pass  the  manufacturer  test 
by  the  number  of  cells  which  include  items  that  have  passed  the  manufacturer  test.  In 
this  example,  we  can  calculate  Mw//  as  follows: 

]  9  0 

•  C  ]  =  2.375  items  for  manufacturer  test 

o.O 

•  [  ]  =  2.500  items  for  temperature  and  humidity  test 

o.O 

•  [  ]  =  0.125  item  for  vibration  test 

o.O 

•  [  ]  =  0.000  item  for  altitude  test 

o.O 

finally  the  expected  number  of  device  for  this  cell  is, 

Mppf}  =  5.000 . 

After  computing  M,lk,  for  each  cell  using  the  same  procedure  above,  we  can  group  these 
cells  to  estimate  Mff00  ,  Mf/O0 ,  M ^  ,  MfJOO  as  follows: 

•  £  Mffil  =  Mff0B 

•  X  M'P'i  =  Mf/00 

•  =  MfpOO 

•  —  M.fpn 

where  k,  /  =  p,  f. 
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In  this  example, 


Table  6.  GRUOPED  DATA  FOR  INITAL  GUESS  CALCULATION 


TERMS 

ESTIMATED  VALUES 

32.500 

22.500 

M/pOO 

23.500 

13.500 

By  combining  the  parameters  n  and  ).'p ,  we  have  four  equations  and  four  unknown 
parameters.  We  can  easily  solve  for  the  parameters  from  the  initial  guesses. 


Table  7.  ESTIMATED  PARAMETERS  FOR  INITIAL  GUESS 


TERMS 

ESTIMATED  VALUES 

n  +  K 

-1.97 

K 

1.08 

K 

1.08 

K 

0.56 

Finally,  we  can  split  the  sum  of  the  n  and  X'p  into  two  parts.  One  possibility  is  to 
divide  the  sum  by  two  assigning  half  to  n  and  half  to  X'p.  Results  using  a  random  initial 
guess  and  the  above  procedure  for  an  initial  guess  are  very  close.  Calculation  of  initial 
guess  parameters  are  done  by  the  FORTRAN  program  INITIAL  in  Appendix  D  and 
the  FORTRAN  program  PARAM  in  Appendix  E. 
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E.  RESULTS 


Reliabilities  of  the  pyrotechnic  device  are  calculated  by  a  FORTRAN  program 
LLMDEP  in  Appendix  F.  They  are  given  in  Table  8.  Results  from  the  log  linear  model 
design  with  dependence  assumption  are  similar  to  the  first  model.  We  can  easily  see 
from  the  previous  table  that  the  ordering  of  the  estimates  of  R  is  counter-intuitive.  We 
expect  (  1  1  2  1  )  to  yield  the  smallest  R  rather  than  an  R  =  0.934  which  is  larger  than 
R  for  about  half  of  the  24  realizations  of  the  failure  vector. 

It  is  clear  from  these  results  that  this  log  linear  model  is  inappropriate  for  modeling 
the  outcomes  of  the  sampling  inspection  plan.  More  realistic  models  would  include 
three-way  and  four-way  interaction  terms.  However,  due  to  the  extreme  censoring  in  the 
data,  we  can  not  estimate  R  for  these  models.  Moreover,  reliabilities  of  several  cases  for 
even  this  model  were  not  calculated  because  of  computational  limitations.  Taking  the 
account  of  dependence  with  log  linear  models  is  clearly  not  a  reasonable  approach  for 
estimating  reliabilities  from  this  data. 
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Table  8.  RELIABILITIES  WITH  LOGLINEAR  MODEL 


FAILURE  VECTOR 

R 

(0000) 

1.0000000 

(1000) 

0.9887247 

(1010) 

0.9693505 

(1001) 

0.9518477 

(1100) 

0.9503552 

(  I  1 I 0) 

0.9457148 

(10  11) 

0.9456805 

(1111) 

0.9399773 

(1020) 

0.9374950 

(1121) 

0.9336501 

(0011) 

0.9255090 

(0110) 

0.9255087 

(1101) 

0.9246104 

(1021) 

0.9175486 

(1120) 

0.9175295 

(0101) 

0.9089977 

(0021) 

0.8978485 

(0120) 

0.8978109 

(0111) 

0.8825021 

(0121) 

0.8584540 
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IV.  WORST  CASE  SCENARIO 
A.  ASSOCIATION  ANALYSIS 

In  the  inspection  sampling  plan,  failures  are  not  assigned  to  an  individual  failure 
mechanism  except  for  the  manufacturer  test.  For  example,  from  the  following  failure 
vector. 


(1  12  1) 

We  know  that  there  is  at  least  one  failure  due  to  the  manufacturing  test.  There  is  one 
failure  from  the  joint  temperature  and  humidity  and  manufacturer  test,  but  we  do  not 
know  which  failure  mechanism  generated  this  failure.  Items  can  fail  due  to  a  manufac¬ 
turing  related  failure  mechanism  or  other  failure  mechanism  or  both.  The  worst  case  is 
to  assume  that  the  cause  of  failure  is  due  to  all  of  the  failure  mechanisms  that  it  was 
exposed  to. 

Let  MT  represents  manufacturer  test,  TAHT  represents  temperature  and  humidity 
test,  VT  represents  vibration  test,  AT  represents  altitude  test. 

Because  all  92  items  are  exposed  to  manufacturer  test,  the  worst  case  of  the  failure 
analysis  of  the  failure  vector  above  is  shown  by  the  following  table. 


Table  9.  FAILURE  ANALYSIS 


MT 

TAHT 

VT 

AT 

#  OF  TESTED  ITEMS 

92 

20 

32 

20 

#  OF  SUCCESFUL 

ITEMS 

87 

19 

30 

19 

This  "  worst  case  "  scenario  should  give  us  lower  bounds  for  the  estimates  R  of  R. 
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If  we  calculate  the  reliability  of  series  system  assuming  failure  mechanisms  in  test 
are  independent  when  in  fact  they  are  associated  but  not  independent,  then  we  underes¬ 
timate  system  reliability.  It  is  reasonable  to  expect  that  there  is  positive  dependence 
between  four  tests.  If  one  item  fails  in  manufacturer  test,  it  is  more  likely  to  fail  in  other 
tests.  One  way  of  modelling  positive  dependance  between  tests  is  to  assume  that  tests 
are  positively  quadrant  dependent. 

Given  random  variables  Tx  ,  ....  Tn.  They  are  said  to  be  positively  quadrant  de- 
pendwm  tPQL>)  1  Kef.  8:  p.  33j,  if 

n 

P(TX  Tn  <  tn)  >  Y\PiTt  <  tt)  (4.1) 

<=i 

for  all  ( tu  t2 ,  .  .  .  ,  r„)  e  R \  An  equivalent  formulation  of  positive  quadrant  depend¬ 
ence  is  7),  T2,  .  .  .  ,  T"  are  PQD  iff 


P{Tt>  tx . T„  >  tH)  >  f]^(  Ti  >  h)  (4.2) 

(=i 

for  all  ( r„  t2,  . . . ,  tn )  e  Rn.  The  proof  that  (4.1)  and  (4.2)  are  equivalent  as  given  in 
(Ref.  8:  pp.  32-33]. 

We  may  take  account  positive  dependence  between  tests  with  the  PQD  assumption. 
Let 


^  = 


1  ,  if  an  item  passes  environment  i 
0  ,  if  an  item  fails  environment  i 


(4.3) 


fori  =  1,2,  3,  4.  Assume  that  7) ,  T2 ,  Tz,  T4  are  positively  quadrant  independent,  i.e. 
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(4.4) 


P(  r,  <  r, ,  t2  <  t2 ,  t3  <  r3 ,  r4  <  f4)  >  f]p(  r,  <  r,) 

i=i 


or  equivalently  that 


P{  Ty  >  /,  ,  T2  >  t2T3  >  t3,  Tt>  tA)  >  ]~]p(  7}  >  /,).  (4.5) 

(=i 

Let  R  be  the  probability  that  an  item  activates  after  exposure  to  all  4  environments. 
Then 


R  =  p{  r,  -  1 ,  t2  =  1 ,  t3  *  1 ,  r4  =  1 ) .  (4.6) 

Using  equations  (4.4)  and  (4.5) 

k  £  P{TX  =  1 )  P(  r2  =  I )  P(  r3  =  1 )  P(  r4  =  1 ) .  (4.7) 

Using  the  notation  from  the  previous  section 

*  ^  g,  &  a  & ,  (4.8) 

where  <2,  is  the  probability  that  an  item  passes  environment  i.  With  the  censored  data, 
we  can't  estimate  g, ,  Q2 ,  g,  and  g,  without  building  a  more  structured  model  for 
T, ,  Tj ,  and  7V  One  alternative  is  get  a  lower  bound  for  estimates  of 
g, ,  gj ,  gj ,  g, .  This  is  worst  case  scenario  . 

B.  CALCULATIONS  WITH  EXAMPLE 
Let 

■  R  be  the  MLE  for  the  reliability  of  the  device. 

■  g,  be  the  MLE  for  the  probability  that  an  item  passes  environment  i. 
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•  R  be  a  lower  bound  of  R  for  reliability  of  the  device. 

■  <2,  be  a  lower  bound  for  the  probability  that  an  item  passes  environment  i. 
Then  according  to  equation  (4.8) 


R  > 


(4.9) 


Further  we  will  construct  estimates  Q,  of  Q,  such  that 


Q,  <  Qi  (4.10) 


and  then  define 


Thus 


A 


i=l 


R  <  R.  (4.12) 
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From  Chapter  II,  the  likelihood  of  observing  Xx  —  x„  X2  =  x2,  X^,  =  x3,  X4  =  x4  is 


x  |"2j/>(r1  =  ij2=i)^ 

x  (1  -  P(TX  =  1,  T2  =  I)  f2~X2  (4.13) 

x  ("3)p<ri“i’r3  =  n* 
x  (i  -p(tx  =  i,r3  =  l))"3-*3 
x  ^jp(tx  =  l,  r4=  i)* 
x  (i  -P(r,  =  i,r4=  i)),,4_^ 


where  x,  is  the  number  of  devices  out  of  n ,  that  activate  after  test  i.  If  we  know  why  the 
device  failed  for  the  tests  i  =  2,  3,  4  which  include  manufacturer  test  along  with  exposure 
to  environment  i  then  our  likelihood  could  be  written  as 
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L(X},  X2,  Xj,  /?2>  P3,  /?4  )  — 


p(t}  =  i)Xi  p(r,  =  o)',i_Xi 


x  p^,  =  i,  r2  =  i)X2P(r,  =  1,  t2  =  0)X21 

x  p (r,  =  o,  r2  =  i)x“  p{t]=o,t2  =  o)X23 


x  p {tx  =  l,  r3  =  i)Xj  p (r,  =  l,  r3  =  o)X31 
x  p (r,  =  o,  r3  =  i)X32  p (r,  =  o,  r3  =  of33 


X 


«4  \ 

>•^4’  **41  *  ***42  >  **43/ 


x  p(r]  =  i,r4=i)^P(r1  =  i,r4  =  o)x-> 
x  p  (r,  =  o,  r4  =  l)^2  p  (r,  =  o,  r4  =  o)^3 


where  jr,,  is  the  number  given  test  i  which  failed  due  to  environment  i  but  passed  manu¬ 
facturing,  jt,2  is  the  number  given  test  i  that  passed  environment  i  but  failed  manufac¬ 
turing  and  xn  is  the  number  given  test  i  that  failed  both  manufacturing  and  environment 
i  for  i  =  2,  3,  4.  Note  that  n,  —  x,  =  xlt  +  xa  +  xn  for  i  =  2,  3,  4. 

From  this  likelihood,  the  MLE's  of  Q,  —  P  {T,=  1)  are  given  by 


Cl 


x,  +  (jc2  4-  .r21)  +  (.r3-f.r31)  +  (x4  +  x41) 
«,  +  n2  +  n3  +  /74 


(4.15) 
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(4.16) 


(4.17) 

(4.18) 

(4.19) 

(4.20) 

(4.21) 

(4.22) 

(4.23) 
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In  Worst  Case  Scenario,  we  assume  that  the  item,  which  fails  in  test  i,  fails  because 
of  both  manufacturing  related  failure  mechanism  and  other  failure  mechanism  which  is 
induced  by  the  i  th  test  environment.  With  this  assumption,  the  reliability  of  of  device 
is  given  as  follows 


R  =  Qi  Qi  Q3  Qa  (4.24) 


where 


Q, 

A 

0.2 

A 

& 

A 

Qa 


X,  +  x2  +  Xj  +  x4 

92 


20 

32 

ii 

20 


Here  is  an  example,  in  this  example  our  failure  vector  is, 

(0  1  2  1  ) 


(4.25) 


^  20  +  19  +  30  +  19 

U'  -  92 

4  =  -^  =  0.9500000 
4=4?-=  0.9375000 
4=4?-=  0.9500000 


0.9456521 


(4-26) 


And  finally  reliability  of  the  device  can  be  estimated  by  using  equation  (4.24) 


R  =  (  0.9456521  )  (  0.9500000  )  (  0.9375000 )  (  0.9500000 )  =  0.8093069  (4.27) 
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C.  RESULTS 


We  compute  MLE's  with  a  FORTRAN  program  MLEB  in  Appendix  G.  These 
are  given  Tables  10  and  11. 

To  get  an  approximate  lower  confidence  bounds  for  R  using  the  worst  case  data, 
we  bootstrap  using  the  procedure  described  in  the  previous  chapter.  The  FORTRAN 
program  RANVEC  in  Appendix  B  is  used  to  generate  the  5000  bootstrap  samples  for 
each  case.  We  then  estimate  R  's  for  5000  failure  vectors  by  the  means  of  MLEB  for 
given  case.  The  next  step  is  to  compute  the  order  statistics  of  R 's  from  /?,  to  Rsooo.  This 
is  done  by  FORTRAN  program  SORT  in  Appendix  C.  Finally,  we  obtain  the  95  % 
lower  confidence  bound  from  this  routine.  Reliabilities  and  95  %  lower  confidence 
bounds  are  tabulated  in  the  following  pages.  Results  of  reliabilities  and  95  %  lower 
confidence  bounds  are  given  in  descending  order. 

These  estimates  and  LCB's  for  R  are  decreasing  with  the  total  number  of  the  failures 
out  of  92  items  tested.  This,  at  least,  is  consistent  with  how  we  believe  pyrotechnic  de¬ 
vices  behave.  It  should  be  noted  that  these  estimates  are  in  fact  conservative  lower 
bounds  for  the  time  MLE's  under  the  very  weak  assumption  of  PQD.  How  conservative 
these  estimates  are  cannot  be  determined  without  more  extensive  data  that  allows  us  to 
estimate  the  degree  of  dependence  between  tests. 
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Table  10.  TEST  PROBABILITIES  (  B  ) 


A 

A 

A  ! 

A 

CASE 

R, 

R , 

R . 

Hum 

1.0000000 

1.0000000 

1.0000000 

1.0000000 

(1000) 

0.9891340 

1.0000000 

1.0000000 

1.0000000 

(0010) 

0.9891340 

1.0000000 

0.9687500 

0.1.00000 

(1010) 

0.9782608 

1.0000000 

0.9687500 

1.0000000 

(0001) 

0.9891304 

1.0000000 

1.0000000 

0.9500000 

KBS3M 

0.9891340 

0.9500000 

1.0000000 

1.0000000 

(1001) 

0.9782608 

1.0000000 

1.0000000 

0.9500000 

mssm i 

0.9500000 

1.0000000 

1.0000000 

1IM1 

0.9782608 

1.0000000 

0.9375000 

1.0000000 

(1020) 

0.9673913 

1.0000000 

0.9375000 

1.0000000 

(001  1  ) 

0.9782608 

1.0000000 

0.9687500 

0.9500000 

(0110) 

0.9782608 

0.9500000 

0.9687500 

1.0000000 

(10  11) 

0.9673913 

1.0000000 

0.9687500 

(1110) 

0.9673913 

0.9500000 

0.9687500 

1.0000000 

(0101) 

1.0000000 

0.9500000 

(1101) 

0.9673913 

0.9500000 

1.0000000 

0.9500000 

0.9673913 

1.0000000 

0.9375000 

0.9500000 

(0120) 

0.9673913 

0.9500000 

0.9375000 

1.0000000 

(1021) 

0.9565217 

1.0000000 

0.9375000 

0.9500000 

(1120) 

0.9565217 

0.9500000 

0.9375000 

1.0000000 

(0111) 

0.9673913 

0.9500000 

0.9687500 

0.9500000 

(Mil) 

0.9565217 

0.9500000 

0.9687500 

0.9500000 

(0121) 

0.9565217 

0.9500000 

0.9375000 

0.9500000 

(1121) 

0.9456521 

0.9500000 

0.9375000 

0.9500000 

Table  11.  RELA1BILITIES  AND  95  %  LOWER  CONFIDENCE 

BOUNDS ( B ) 


FAILURE 

VECTOR 

A 

R 

FAILURE 

VECTOR 

95  %  LCB 

(0  0  0  0  ) 

1.0000000 

1.0000000 

(  10  0  0) 

0.9891304 

(  1  0  0  0  ) 

0.9891304 

(0010) 

0.9582200 

(0010) 

0.8766983 

(1010) 

0.9476901 

(1010) 

0.8766983 

(  0  0  0  1  ) 

0.9396738 

(  0  0  0  1  ) 

0.8222825 

(0100) 

0.9396738 

— na— 1 

0.8222825 

(10  0  1) 

0.9293478 

m 

0.8222825 

(1  10  0) 

0.9293478 

(1100) 

0.8222825 

(  0  0  2  0  ) 

0.9171195 

■unit 

0.8179346 

(  1  0  2  0  ) 

0.9069293 

■n mm 

0.8084239 

(0011) 

0.9003056 

M!M— 

0.7712974 

(0110) 

0.9003056 

(0110) 

0.7712974 

(1011) 

0.8903022 

(1011) 

0.7712974 

(1110) 

0.8903022 

(0101) 

0.7712974 

(0101) 

0.8828803 

(0101) 

0.7548369 

(1101) 

0.8730705 

(1101) 

0.7548369 

(0120) 

0.8615828 

(0120) 

0.7218070 

(0021) 

0.8615828 

(0021) 

0.7200747 

(1021) 

0.8519021 

(1021) 

0.7200747 

(1120) 

0.8519021 

(1120) 

0.7200747 

(0111) 

0.8457870 

(0111) 

0.6929346 

(1111) 

0.8362838 

(1111) 

0.6927614 

(0121) 

0.8093069 

(0121) 

0.6508338 

(112  1) 

0.8001103 

(1121) 

0.6505433 
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V.  BONUS  SYSTEM  APPROACHES 

A.  BACKGROUND 

Quality  is  described  as  especially  high  degree  of  goodness  or  worth  "  [Ref.  9 
:  p.  685].  In  industry,  a  quality  product  is  one  that  fulfills  customer  expectations. 
There  are  two  general  aspects  of  quality: 

■  Quality  of  Design 

■  Quality  of  Conformance. 

All  goods  and  services  are  produced  in  various  grades  or  levels  of  quality.  These 
variations  in  grades  or  levels  of  quality  are  intentional;  therefore  the  appropriate  tech¬ 
nical  term  is  quality  of  design.  The  quality  of  conformance  is  how  well  the  product 
conforms  to  the  specifications  and  tolerances  required  by  the  design.  Quality  of  con¬ 
formance  is  influenced  by  a  number  of  the  following  factors: 

■  the  choice  of  manufacturing  process 

■  the  training  and  supervision  of  the  workforce 

■  type  of  quality  assurance  system  (  process  controls,  tests,  inspections,  etc.  )  used 

■  the  extent  to  which  these  quality-assurance  procedures  are  followed 

■  the  motivation  of  workforce  to  achieve  quality. 

Quality  Control  is  the  engineering  and  management  activity  by  which  we  measure 
the  quality  charactenstics  of  a  product,  comparing  them  with  specifications  or  require¬ 
ments  and  taking  appropriate  remedial  action  whenever  there  is  a  difference  between  the 
actual  performance  and  the  standard  [Ref.  10  :  pp.  1-3). 

As  contracts  are  now  written,  contractors  need  only  to  satisfy  the  requirments  of  the 
sampling  inspection  plan  for  lot  acceptance.  Contractors  have  no  incentive  to  improve 
the  quality  of  the  items  they  provide,  although  they  are  in  a  position  to  do  so.  As 
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mentioned  before,  the  quality  of  conformance  is  influenced  by  the  motivation.  Therefore 
we  can  motivate  manufacturers  by  giving  a  bonus  for  improved  quality.  To  improve 
quality,  The  Naval  Weapons  Support  Center  has  decided  to  implement  a  bonus  system 
for  pyrotechnic  devices. 

B.  BONUS  PLANS 

In  this  chapter,  we  design  a  bonus  system  to  improve  the  quality  of  pyrotechnic 
devices.  A  good  bonus  system  encourages  reliability  growth,  because  firms  try  to  reach 
a  high  quality  to  get  a  bonus.  An  effective  bonus  system  must  detect  small  differences 
among  the  offered  lots.  In  the  previous  chapters,  we  estimated  the  reliability  of  a 
pyrotechnic  device  in  several  ways.  We  assumed  independence  in  the  first  model,  and 
we  assumed  dependence  between  tests  in  the  last  two  models.  The  estimated  reliabilities 
in  the  first  two  models  are  close  to  each  other,  but  they  exhibit  different  structures  in 
order.  One  can  easily  see  that  the  failure  vector  of  the  worst  case  (112  1)  has  bigger 
reliability  value  than  ten  of  the  possible  cases,  from  Table  5  and  Table  15.  In  addition, 
the  order  of  95  %  LCB's  of  cases  does  not  match  to  the  order  of  MLE's  in  independent 
models.  But  in  the  worst  case  scenario,  which  assumes  dependence  between  tests,  we 
get  the  same  order  for  both  95  %  LCB's  and  MLE's.  Thus,  we  will  use  the  worst  case 
scenario  model  in  further  calculations  for  the  bonus  system.  With  the  assumptions 
above,  ure  can  apply  three  sampling  plans  for  giving  bonuses  to  the  manufacturers. 
These  are: 

*  Single  Sampling  Bonus  System 

■  Double  Sampling  Bonus  System 

■  Multi-Sampling  Bonus  System. 

Manufacturers  have  to  meet  pyrotechnic  device  acceptance  criteria  first,  before 
having  a  chance  to  get  a  bonus. 
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1.  Single  Sampling  Bonus  System 

In  this  bonus  system,  we  will  decide  two  things  at  the  end  of  inspection.  First 
we  will  decide  whether  the  ofTered  lot  is  acceptable  or  not,  and  then  for  lots  which  are 
accepted  we  will  decide  whether  to  give  a  bonus  to  the  manufacturer. 

The  Single  sampling  bonus  system  is  designated  by  three  numbers.  These  are, 

w,  ,  LCBFB  ,  LCB 


where 

■  n,  is  the  sample  size  for  inspection  (  here  «,  =  92  ). 

■  LCBFB  is  the  cut  off  value  applied  to  the  lower  confidence  bound  for 
awarding  a  bonus. 

■  LCB  means  estimated  lower  confidence  bound  after  inspection.  (If  LCB  is 
greater  than  or  equal  to  LCBFB,  then  a  bonus  is  awarded.) 

The  following  algorithm  shows  us  how  single  sampling  bonus  system  works  for 
pyrotechnic  devices. 


STEP  #  1 
STEP  #  2 
STEP  #  3 
STEP  #  4 
STEP  #  5 

STEP  #  6 
STEP  #  7 


DETERMINE  CUT  OFF  VALUE  FOR  AWARDING  THE  BONUS. 
TAKE  A  SAMPLE  SIZE  OF  92  FROM  OFFERED  LOT. 

APPLY  MANUFACTURER  AND  THREE  ENVIRONMENT  TESTS. 
COMPARE  RESULTS  OF  TESTS  WITH  ACCEPTANCE  CRITERIA. 
ESTIMATE  ITS  LOWER  CONFIDENCE  BOUND  (  LCB  ), 

IF  THEY  MET  ACCEPTANCE  CRITERIA. 

COMPARE  LCB  OF  LOT  WITH  LCBFB 

IF  THE  LCB  IS  GREATER  THAN  OR  EQUAL  TO  LCBFB, 

GIVE  BONUS  TO  THE  FIRM. 
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According  to  the  algorithm  described  above,  the  following  are  possible  events 
for  the  firms. 

■  Firm  may  not  satisfy  our  acceptance  criteria.  This  means  that  firm  gets  a 
failure  vector  worse  than  (112  1). 

■  Firm  satisfies  acceptance  criteria,  but  its  LCB  may  be  less  than  LCBFB.  This 
means  that  firm  does  not  get  bonus,  but  the  lot  is  accepted. 

■  Firm  satisfies  acceptance  criteria,  and  its  LCB  may  be  greater  than  or  equal 
to  LCBFB.  This  means  that  firm  gets  the  bonus. 

2.  Double  Sampling  Bonus  System 

A  double  sampling  plan  has  an  advantage  over  a  single  sampling  plan.  Because 
a  double  sampling  plan  involves  a  larger  sample  size,  it  reduces  the  chance  that  a  man¬ 
ufacturer  who  deserves  a  bonus  will  not  get  one.  Double  Sampling  Bonus  System  per¬ 
mits  the  taking  of  two  samples  on  which  to  make  a  decision  [Ref.  11:  pp.  184-185). 

In  this  system,  we  have  two  inspection  stages.  If  the  firm  does  not  get  a  bonus 
after  the  first  inspection,  then  the  firm  is  given  a  second  chance  with  a  second  inspection. 
A  double  sampling  bonus  system  is  designated  by  five  numbers. 

«,  ,  n2 ,  LCBFB ,  LCBl  ,  LCB2 


where 

■  n,  is  the  sample  size  for  first  inspection,  (  here  «,  =  92  ) 

■  n2  is  the  sample  size  for  second  inspection,  (  here  n2  =  92  ) 

■  LCBFB  is  the  cut  off  value  applied  to  lower  confidence  bounds  for  awarding 

bonuses, 

*  LCBl  is  the  estimated  lower  confidence  bound  after  first  inspection, 

■  LCB2  is  the  updated  estimate  of  lower  confidence  bound  after  second  in¬ 
spection. 
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In  the  double  sampling  bonus  system,  if  the  firm  does  not  get  a  bonus  after  the 
first  inspection  but  does  meet  the  acceptance  criteria  for  the  first  sample,  then  a  second 
sample  is  taken.  We  already  have  LCB's  of  possible  cases  aftci  the  first  inspection  in 
Table  1 1.  After  the  second  inspection,  we  calculate  the  LCB  using  an  aggregated  failure 
vector  which  includes  failures  from  both  samples.  When  we  compute  LCB's  using  ag¬ 
gregated  failure  vector  after  the  second  inspection,  we  have  a  total  sample  size  of  1S4; 
there  are  135  different  failure  vectors  for  which  the  lot  meets  the  acceptance  criteria  for 
both  samples.  After  tabulating  these  possible  135  cases  and  the  estimates  of  R,  we  used 
the  bootstrap  procedure  to  find  LCB's  for  each  case.  We  created  5000  random  failure 
vectors  for  each  of  them  by  using  case  success  probabilities  in  tests.  After  this,  we  esti¬ 
mated  YlLE's  of  5000  failure  vectors  for  each  possible  case.  Finally,  we  estimated  95 
%  LCB  of  each  case.  These  calculations  were  done  by  using  programs  RANVEC, 
MLEB,  SORT  in  appendix  B,  G  and  C  respectively.  The  results  are  tabulated  at  Ap¬ 
pendix  I.  The  following  algorithm  shows  us  how  double  sampling  bonus  system  works 
for  pyrotechnic  device. 


STEP  #  1  : 

STEP  #  2  : 

STEP  it  3  : 
STEP  it  4  : 
STEP  it  5  : 

STEP  it  6  : 
STEP  it  7  : 

STEP  it  8  : 


DETERMINE  LCBFB . 

TAKE  A  SAMPLE  SIZE  OF  92  FROM  OFFERED  LOT. 

APPLY  MANUFACTURER  AND  THREE  ENVIRONMENT  TESTS. 
COMPARE  RESULTS  OF  TESTS  WITH  ACCEPTANCE  CRITERIA. 
ESTIMATE  ITS  LOWER  CONFIDENCE  BOUND  (  LCB1  ), 

IF  THEY  MET  ACCEPTANCE  CRITERIA. 

COMPARE  LCB1  OF  LOT  WITH  LCBFB 

IF  THE  LCB1  IS  GREATER  THAN  OR  EQUAL  TO  LCBFB, 

GIVE  BONUS  TO  THE  FIRM. 

IF  THE  LCB1  IS  LESS  THAN  LCBFB, 
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GIVE  A  SECOND  CHANCE  TO  THE  FIRM  FOR  BONUS. 


STEP  #  9 

STEP  #  10 
STEP  #  11 

STEP  #  12 

STEP  #  13 
STEP  //  14 


TAKE  A  NEW  SAMPLE  SIZE  OF  92  FROM  OFFERED  LOT. 
APPLY  MANUFACTURER  AND  THREE  ENVIRONMENT  TESTS. 
AFTER  GETTING  THE  NEW  FAILURE  VECTOR, 

ADD  THIS  ONE  TO  THE  FIRST  FAILURE  VECTOR. 
ESTIMATE  ITS  LOWER  CONFIDENCE  BOUND  (  LCB2  ), 
WITH  AGGREGATED  FAILURE  VECTOR. 

COMPARE  LCB2  OF  LOT  WITH  LCBFB 

IF  THE  LCB2  IS  GREATER  THAN  OR  EQUAL  TO  LCBFB, 

GIVE  BONUS  TO  THE  FIRM. 


According  to  the  algorithm  described  above,  the  following  are  possible  events 
for  the  firms. 

■  The  firm  may  not  satisfy  our  acceptance  criteria  after  first  inspection.  This 
means  that  firm  gets  a  failure  vector  worse  than  (112  1).  The  lot  is  not  accepted  and 
second  sample  is  not  taken. 

■  The  firm  satisfies  acceptance  criteria,  and  its  LCB1  may  be  greater  than  or 
equal  to  LCBFB.  This  means  that  it  gets  bonus  after  first  inspection,  and  a  second 
sample  is  not  needed. 

■  The  firm  satisfies  acceptance  criteria,  but  its  LCB1  may  be  less  than  LCBFB. 
This  means  that  it  does  not  get  the  bonus  after  first  inspection,  but  still  has  a  chance  to 
get  a  bonus  if  it  submits  a  second  sample. 

Of  those  which  submit  a  second  sample 

■  LCB2  may  be  less  than  LCBFB.  It  means  that  firm  does  not  get  the  bonus, 
but  the  lot  is  accepted. 
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■  The  aggregated  failure  vector  meets  acceptance  criteria.  But  LCB2  may  be 
grater  than  or  equal  to  LCBFB.  It  means  that  firm  gets  bonus. 

3.  Multi-Sampling  Bonus  System 

In  this  system  approach,  we  may  take  a  sample  of  size  9?  from  offered  lot  three 
times,  four  times  or  more.  But  if  we  use  the  three  sampling  bonus  system,  after  the  third 
inspection,  we  may  get  (  4  X  4  X  7  x  4  =  44S  )  448  different  failure  vectors.  Multi- 
Sampling  Bonus  System  is  going  to  be  more  computationally  intensive  with  respect  to 
double  sampling  bonus  system.  This  can  be  done  later  using  with  the  same  reasoning 
in  double  sampling  bonus  system. 

C.  EXAMPLES 

1.  Single  Sampling  Bonus  System 

Let  us  assume  that  we  decided  that  lower  confidence  bound  for  bonus  will  be 
0.8500.  We  are  going  to  use  this  LCBFB  in  these  examples. 

■  Case  #  1 

At  the  end  of  inspection,  firm  has  following  failure  vector 

(  2  0  0  0  ). 

The  firm  did  not  meet  the  acceptance  criteria.  We  immediately  reject  the  lot  and  no 
bonus  is  given. 

■  Case  #  2 

At  the  end  of  inspection,  firm  has  following  failure  vector 

(10  0  1  ). 

The  firm  meets  the  acceptance  criteria.  We  estimate  its  95  %  lower  confidence  bound 
to  be 
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LCB  =  0.8222825  . 


Because  LCB  <  LCBFB  ,  we  do  not  give  a  bonus  to  the  firm. 

■  Case  #  3 

At  the  end  of  inspection,  firm  has  following  failure  vector 

(001  0  ). 

The  firm  meets  the  acceptance  criteria.  We  estimate  its  95  %  lower  confidence  bound 
to  be 


LCB  =  0.87669S3  . 

Because  LCB  >  LCBFB  ,  we  give  a  bonus  to  the  firm. 

2.  Double  Sampling  Bonus  System 

•  Case  #  1 

At  the  end  of  first  inspection,  firm  has  following  failure  vector. 

(10  0  1  ). 

The  firm  meets  the  acceptance  criteria.  We  estimate  its  95  %  lower  confidence  bound 
to  be 


LCB  1  =  0.8222825  . 

Because  LCB  1  <  LCBFB  ,  we  give  a  second  chance  to  the  firm.  Here  is  the  result  of 
second  inspection 


(  0  3  2  1  ). 

The  aggregated  failure  vector  will  be 
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(  1  3  2  2  ) . 


It  is  clear  that  LCB2  <  LCB1  for  this  data,  thus  we  do  not  give  a  bonus  to  the  firm. 

■  Case  #  2 

At  the  end  of  first  inspection,  firm  has  following  failure  vector 

(100  1  ). 

The  firm  meets  the  acceptance  criteria.  We  estimate  its  95  %  lower  confidence  bound 
to  be 

LC2?1  =  0.8222825  . 

Here,  LCB  <  LCBFB  thus  we  will  give  a  second  chance  to  the  firm.  Here  is  the  result 
of  second  inspection 

(001  1  ). 

The  aggregated  failure  vector  will  be 

(10  1  2  ) , 

and  we  estimate  its  95  %  lower  confidence  bound  to  be 

LCB\  =  0.8369564  . 

Because  LCB  <  LCBFB  ,  we  do  not  give  a  bonus  to  the  firm. 

■  Case  #  3 

At  the  end  of  first  inspection,  firm  has  following  failure  vector 

(10  0  1  ). 
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The  firm  meets  the  acceptance  criteria.  We  estimate  its  95  %  lower  confidence  bound 
to  be 

LCBl  =  0.8222825  . 

Because  LCBl  <  LCBFB  ,  we  will  give  a  second  chance  to  the  firm.  Here  is  the  result 
of  second  inspection 

(  0  0  0  1  ). 

The  aggregated  failure  vector  will  be 

(  1  0  0  2  ) , 

and  we  estimate  its  95  %  lower  confidence  bound  to  be 

LCBl  =  0.8706521  . 

Because  LCBl  >  LCBFB  ,  we  will  give  a  bonus  to  the  firm. 
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VI.  SIMULATION  RESULTS  OF  BONUS  SYSTEMS 

A.  BACK  GROUND 

We  simulated  the  bonus  system  to  get  an  idea  of  how  well  the  bonus  system  pro¬ 
posed  in  the  previous  chapter  works.  This  provides  the  user  with  a  means  of  setting  the 
cut  off  criteria  for  awarding  a  bonus.  To  generate  random  failure  vectors,  we  need  to 
know  the  probabilities  of  being  successful  in  each  test  for  the  firm.  First,  we  assume  that 
the  firm  has  an  equal  probability  of  being  successful  in  each  test,  with  the  following 
values  of  P,  =  I\  ~  Pi  =  P4 

•  0.9375 

•  0.9500 

•  0.9750 

•  0.9900 

•  0.9950. 

In  a  second  set  of  simulations,  we  assume  that  the  manufacturer  test  has  a  bigger 
probability  of  being  successful  than  the  other  environment  tests.  The  following  proba¬ 
bilities  of  being  successful  in  manufacturer  test  were  used. 

•  0.9990 

•  0.9950 

•  0.9750 

•  0.9500. 
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In  these  cases,  we  used  equal  probabilities  of  being  successful  in  other  environment  tests 
with  values  between  0.9375  and  manufacturer  test  probability  for  that  given  case. 
Finally,  we  assume  that  the  manufacturer  test  should  have  smaller  probability  of  being 
successful.  When  we  study  these  cases  with  this  assumption,  we  used  probabilities  of 
being  successful  in  other  environment  tests  between  manufacturer  test  probability  and 
0.9990  for  that  given  case.  We  assume  that  minimum  value  of  the  manufacturer  test 
probability  will  be  0.9500  because  of  worst  case  in  last  two  assumptions. 

We  generated  2000  random  failure  vectors  for  each  possible  combination  of  proba¬ 
bilities  by  using  program  RANVEC  in  Appendix  B.  They  were  used  for  1000  repli¬ 
cations  of  each  bonus  system,  because  the  Double  Sampling  Bonus  System  (  DSBS  ) 
potentially  uses  two  failure  vectors  per  replication. 

After  getting  the  failure  vectors,  the  next  step  is  to  decide  lower  confidence  bound 
for  giving  bonus  (  LCBFB  ).  We  chose  0.800,  0.825,  0.850,  0.875,  0.900,  0.950,  0.999 
as  LCBFB  during  our  simulations.  We  simulated  bonus  systems  by  using  program  BO¬ 
NUS  in  Appendix  G.  Program  BONUS  counts  how  many  times  firm  gets  bonus  during 
1000  replications.  And  then  it  calculates  the  bonus  percentage  dividing  counted  number 
by  1000.  For  each  scenario  this  bonus  percent  is  an  estimate  of  the  probability  of  getting 
a  bonus. 

B.  INITIAL  COMPARISON  OF  SYSTEMS 

In  this  part,  with  equal  probabilities  of  being  successful  in  tests,  we  tried  to  see  the 
dilfercnce  between  the  Single  and  the  Double  Sampling  Bonus  Systems.  For  this  reason, 
we  simulated  Single  Sampling  Bonus  System.  Results  were  tabulated  and  plotted  in  next 
pages. 
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Table  12.  SINGLE  SAMPLING  BONUS  SYSTEM  (EQUAL  PROBABILITIES) 


PROB.'S 

CHOSEN  LOWER  CONFIDENCE  BOUNDS  FOR  BONUS 

0.800 

0.825 

0.850 

0.875 

0.900 

0.950 

0.999 

0.9200 

0.025 

0.007 

0.007 

0.007 

0.002 

0.002 

0.002 

0.9375 

0.060 

0.023 

0.023 

0.023 

0.010 

0.010 

0.007 

0.9450 

0.092 

n 

0.044 

0.044 

0.018 

0.018 

0.010 

0.9500 

0.117 

0.053 

0.053 

0.053 

0.022 

0.022 

0.011 

0.9600 

0.222 

0.112 

0.045 

B 

0.028 

0.9700 

0.380 

m 

IB 

0.198 

0.107 

0.107 

0.067 

0.9750 

U.482 

m 

H| 

■ 

0.163 

0.067 

0.9800 

KS3 

m 

0.381 

H| 

B 

0.237 

0.106 

0.9850 

m 

m 

m 

B 

B 

B 

0.166 

0.9900 

m 

m 

Q 

B 

B 

B 

0.274 

mm 

b 

m 

0.805 

m 

0.502 

B 

B 

It  is  obvious  that  there  is  no  difference  between  some  lower  confidence  bounds  for  bonus 
from  the  table  above.  The  firm  gets  the  same  bonus  percentage  when  we  use  0.825, 
0.850,  0.875  as  lower  confidence  bound  for  bonus.  The  same  thing  occurs  when  we  use 
0.900  and  0.950  as  LCBFB.  For  this  reason,  we  are  going  to  see  four  curves  in  Figure 
1.  Figure  2  shows  Double  Sampling  Bonus  System  with  different  LCBFB's. 
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Figure  1.  Single  Sampling  Bonus  System  With  Different  LCBFB's 


DOUBLE  SAMPLING  BONUS  SYSTEM 
DIFFERENT  LCB'S  FOR  BONUS  WfTH  EQUAL  TEST  PROBABILITIES 


0-12  OJK  O.SS  0.96 


PROBABILITIES 

Figure  2.  Double  Sampling  Bonus  System  With  Different  LCBFB's 
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As  you  sec  in  figure  2,  we  have  six  different  curves  for  Double  Sampling  Bonus  System. 
Because  in  this  system,  we  can  see  from  the  table  that  LCBFB's  0.900  and  0.950  have 
approximately  the  same  bonus  percentage,  f  or  this  reason  we  did  not  plot  for  0.950. 


Table  13.  DOUBLE  SAMPLING  BONUS  SYSTEM  (EQUAL  PROBABILITIES) 


PROB.'S 

CHOSEN  LOWER  CONFIDENCE  BOUNDS  FOR  BONUS 

m 

B 

B 

0.875 

0.900 

B 

0.999 

0.9200 

0.026 

0.008 

B 

0.007 

B 

0.002 

0.002 

0.9375 

0.064 

0.030 

0.025 

0.023 

0.010 

B 

0.007 

■ 

0.9450 

0.108 

0.059 

m 

0.045 

0.020 

0.018 

0.010 

0.9500 

m 

B 

0.057 

B 

n 

0.011 

0.9600 

0.277 

0.172 

129 

H 

0.045 

0.028 

0.9700 

B 

B 

0.293 

0.236 

B 

0.067 

0.9750 

0.624 

B 

m 

B 

0.212 

m 

0.106 

0.9800 

0.752 

m 

B 

m 

B 

0.237 

0.9850 

B 

m 

B 

m 

B 

m 

B 

■ ■ 

B 

B 

B 

m 

B 

m 

0.407 

0.9950 

B 

B 

B 

B 

B 

Figure  2  shows  DSBS  to  be  more  sensitive  in  the  sense  that  a  higher  percentage  of 
the  firms  were  awarded  a  bonus.  For  this  reason,  we  decided  to  implement  DSBS  in  all 
simulations. 
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C.  SIMULATION  RESULTS  WITH  DIFFERENT  LCB'S  FOR  BONUS 

In  this  section,  results  for  each  chosen  lower  confidence  bound  for  bonus  will  be 

presented  as  follows.  We  used  values  in  Table  13  to  draw  plots  with  equal  probabilities 

in  each  test.  Polynomial  approximation  was  used  in  curve  fitting. 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  -  0.B00 


Figure  3.  Double  Sampling  Bonus  System  With  LCBFB  =  0.800 


Table  14.  DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.800 


LOW  ER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.800 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

mm 

mm 

mm 

0.9500 

mm 

mm 

mm 

mm 

0.9750 

mm 

mm 

mm 

0.9850 

0.861 

0.9900 

mm 

0.9950 

0.987 
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In  Figure  3,  when  the  probabilities  of  being  successful  in  each  test  increase,  then  the 
bonus  percent  increases.  For  example,  when  the  test  probabilities  is  equal  to  0.9200, 
then  the  bonus  percent  is  0.026.  If  the  firm  increases  its  probabilities  of  being  succesful 
in  each  test  to  0.9900,  then  the  bonus  percent  becomes  0.948.  The  result  of  the  double 
sampling  bonus  system  with  different  lower  confidence  bounds  are  tabulated  and  plotted 
in  Appendix  J. 

To  see  whether  the  probability  of  passing  the  manufacturer  test  effects  the  bonus 
percentage  differently  than  the  other  probabilities,  we  assumed  that  the  environment 
tests  would  have  equal  probabilities.  But  the  firm  will  have  a  different  probability  of 
being  successful  in  the  manufacturer  test.  In  the  following  table,  the  first  row  represents 
probability  of  being  successful  in  the  environment  tests  and  the  first  column  represents 
probability  of  being  successful  in  the  manufacturer  test.  The  intersection  of  rows  and 
columns  gives  us  the  bonus  percentage  of  a  firm  with  given  probabilities  of  being  suc¬ 
cessful  in  the  tests.  These  procedures  were  done  for  each  LCBFB  value  separately.  The 
results  are  tabulated  and  plotted  in  Appendix  K.  Bonus  percentages  for  LCBFB  = 
0.8000  are  plotted  on  the  following  page.  It  is  clear  that  the  bonus  percentage  of  the 
firm  will  be  high  if  the  firm  has  big  probabilities  of  being  succesful  in  both  manufacturer 
test  and  other  joint  environment  tests. 
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D.  BONUS  PERCENTAGE  (  BPRCT  )  FORMULATION 


In  this  section  we  try  to  approximate  probability  of  getting  a  bonus  as  a  function 
of  lower  confidence  bound  for  bonus  (  LCBFB  ),  probability  of  passing  manufacturer 
test  (  P,  )  ,  and  probability  of  passing  environment  tests  (  P2  ).  To  do  this,  we  use  re¬ 
gression  analysis  with  GRAFSTAT  using  the  bonus  percentages  which  are  the  simulated 
values  of  the  probability  of  getting  a  bonus.  If  a  reasonable  relationship  between  the 
bonus  percentage  and  LCBFB  and  the  probabilities  P,  and  P2  can  be  found;  it  can  be 
used  to  set  LCBFB  without  resorting  to  simulation. 

We  want  to  formulate  bonus  percentages  as  a  function  of  the  following: 

■  Lower  confidence  bound  for  bonus  (  LCBFB  ) 

■  Probability  of  passing  manufacturer  test  P, 

■  Probability  of  passing  environment  tests  P2. 

After  polynomial  approximation,  we  can  see  that  plotted  graphs  (Figures  3-10)  look 
like  logistic  growth  curves  [Ref.  12:  p.  383].  Then  we  can  formulate  bonus  probabilities 
(BP)  using  the  logistic  growth  function. 


BP 


1 


1  -f  e 


A 


(7.4) 


where 

A  =  p0  +  /?,  P,  +  p2  P]  +  ft  P2  +  P\  +  05  P,  P2  +  f}6  LCBFB  (7.5) 


We  can  make  a  transformation  described  as  below. 


1  +  e' 


BP 


a  1  -  BP 
e  = 


BP 


(7.6) 


(7.7) 
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.  ,  r  i  -  bp  ' 

A  08  L  bp 


(7.8) 


Now  we  have  a  linear  equation  as  a  function  of  LCBFB,  Px  and  P2.  We  may  use 
simulation  results  and  approximate  parameters  doing  a  linear  regression. 

■  Let  Y  be  the  (  n  x  1)  column  vector  of  observations  on  dependent  variable. 

■  Let  X  be  the  (  n  x  p’  )  matrix  consisting  of  a  column  of  ones,  which  is  labeled  1, 
followed  by  p  column  vectors  of  observations  on  independent  variables. 

■  Let  p  be  the  (  p‘  x  1  )  vector  of  parameters  to  be  estimated. 

■  Let  £  be  the  (  n  x  1  )  vector  of  random  errors.  Then 

Y=Xp  +  c  (7.9) 

We  can  obtain  observations  on  dependent  variable  as  below, 

„  ,  T  1  -  BPRCT  1 

Y  -  loH  bprcT  J  (7  ,0) 

where  BPRCT  is  the  bonus  percentage  obtained  from  the  simulation.  As  an  example 
we  have  the  following  information  from  simulations. 

•  LCBFB  =  0.875 

•  P,  =  0.9500 
.  P\  =  0.9025 

•  P2  =  0.9375 
.  P\  =  0.8789 

.  />,  P2  =  0.8910 
.  BPRCT  =  0.029 
.  Y  =  3.5110. 

We  can  write  following  equation. 
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3.5110  =  P0  +  Pi  0.95  +  p2  0.9025  +  /?3  0.9375  +  /?4  0.8789  +  05  0.8910  +  pt  0.875  +  c. 


We  deliberately  chose  30  random  results  from  simulations.  We  formulate  them  as 
above.  We  can  estimate  unknown  parameters  on  these  equations  doing  linear  re¬ 
gression.  30  equations  can  be  written  with  matrix  notation  as  follows: 


>'d) 

>’(2) 

r(3) 


>'(30) 


(7.11) 


And  X  matrix  is  designed  as  follows: 


1  Pl(l)  p]{  1) 

Pi( » 

^(U 

Pi  Pad) 

1  Pp  2)  P]{  2) 

Pi(2) 

P\{  2) 

^(2) 

X  = 

1  />,(  3)  Pi(  3) 

PjO) 

.  .  ^ 

Z5,  P2(3) 

1  PP 30)  Pi  (30) 

P2(  30) 

P22(30) 

Z5,  Z>2(30) 

And  parameter  vector  will  be  as  below. 


P  = 


Po 

Pi 

Pi 

Pi 

P<i 

Ps 

Pe 


(7.13) 
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We  estimate  parameters  with  linear  regression  using  GRAFSTAT  packages.  And 
then  we  formulate  the  bonus  percentage  using  the  estimates  of  parameters.  The  esti¬ 
mated  bonus  probabilities  are  given  by 


BP  = 


1  +  e 


(7.14) 


l'he  following  estimates  of  the  parameters  were  obtained  from  the  regression. 


Table  16.  COEFFICIENT  V  ECTORS  OF  REGRESSION  ANALYSIS 


PARAMETERS 

ESTIMATED  VALUE 

Po 

-1042.3 

p> 

695.11 

p> 

-86.661 

p* 

1560 

Ps 

-569.81 

p* 

-570.62 

P 7 

12.266 

The  standard  error  was  0.45471  after  the  linear  regression.  We  can  use  these  esti¬ 
mates  of  parameters  in  equation  (7.5).  And  we  can  estimate  the  bonus  probability  of  the 
firm  when  the  probability  of  passing  tests  and  lower  confidence  for  bonus  are  known. 
Some  examples  are  calculated  with  equation  (7.5).  Results  of  the  calculations  and 
comparison  with  simulation  are  summarized  in  following  table. 
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Table  17.  BONUS  PERCENTAGES  WITH  REGRESSION  ANALYSIS 


LCBFB 

/», 

SIMU¬ 

LATION 

EST.  VALUE 

0.800 

0.9750 

0.9750 

0.624 

0.564 

0.800 

0.9990 

0.9375 

0.106 

0.062 

0.825 

0.9750 

0.9950 

0.882 

0.911 

0.825 

0.8500 

0.9750 

0.064 

0.077 

0.850 

0.9500 

0  9900 

0.612 

0.615 

0.850 

0.9900 

0.9750 

0.476 

0.529 

0.875 

0.9000 

0.9750 

0.1 14 

0.080 

0.875 

0.9700 

0.9700 

0.236 

0.211 

0.900 

0.9950 

0.9990 

0.985 

0.945 

0.900 

0.8750 

■SB 

0.158 

0.096 

0.950 

eg 

mm 

0.018 

0.010 

0.950 

EE 

mm 

B 

0.790 

0.999 

B 

0.9850 

0.274 

0.999 

mmm 

B 

0.004 

0.002 

i 
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VII.  CONCLUSIONS  AND  RECOMMENDATIONS 

Wc  estimated  the  reliability  of  pyrotechnic  device  from  the  sampling  plan  using  the 
following  models: 

*  Maximum  likelihood  assuming  independence  of  pyrotechnic  device  activation  in 
different  environments; 

•  Log  linear  model  incorparating  some  dependence  between  pyrotechnic  device  ac¬ 
tivation  in  different  environments; 

■  Worst  case  scenario  which  gives  a  lower  bound  for  the  estimated  reliability  for  the 
general  model  where  no  assumptions  are  made  about  the  form  of  dependence  between 
pyrotechnic  device  activation  in  different  environments. 

Using  these  models  and  based  on  sampling  plan  data,  estimates  of  overall  reliability 
along  with  95  %  lower  confidence  bound  are  obtained.  We  computed  the  lower  confi¬ 
dence  bound  for  each  possible  case  by  bootstrapping. 

Results  from  the  first  model  are  not  consistent  with  the  way  pyrotechnic  devices 
operate.  In  particular,  the  estimated  reliabilities  are  not  ordered  as  we  expected  them 
to  be.  Intiutively,  we  expect  samples  with  fewer  failures  to  give  smaller  R  values  and 
corresponding  LCB's  than  samples  with  more  failures.  This  is  not  the  case  for  the  first 
model.  For  example,  the  failure  vector  that  has  the  maximum  number  of  failures  for 
each  test  (  a  total  of  5  failures  )  is  in  the  middle  of  the  order  with  respect  to  R  and  the 
95  %  lower  confidence  bound.  The  failure  vector  (0101)  with  a  total  of  2  failures 
has  a  much  lower  R  and  LCB.  The  discrepency  between  the  results  of  this  model  and 
what  we  expected  to  see  are  probably  due  to  the  fact  that  there  is  a  dependence  between 
the  events  which  a  device  activates  under  different  environments. 
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The  results  of  the  log  linear  model  design  with  two-way  dependence  assumptions  are 
similar  to  the  first  model.  Two-way  interaction  terms  were  used  in  this  model.  More¬ 
over,  reliabilities  of  several  cases  were  not  calculated  because  of  computational  limita¬ 
tions.  Due  to  the  extensive  censoring  in  the  sample  data,  it  does  not  appear  to  be 
possible  to  estimate  the  reliability  based  on  models  which  incorparate  dependence. 
Thus,  we  turn  to  finding  lower  bounds  for  the  estimated  reliability  based  on  models  with 
dependence. 

Finally,  the  worst  case  scenario  model  gives  the  most  reasonable  results  for  this 
problem.  Both  the  lower  bounds  for  estimates  of  reliabilities  and  95  %  lower  confidence 
bounds  for  these  lower  bounds  are  ordered  according  to  the  total  number  of  failures. 
Thus,  the  results  from  the  worst  case  scenario  were  used  to  implement  the  bonus  system 
for  pyrotechnic  devices. 

After  getting  the  estimate  of  reliability  and  95  %  lower  confidence  bound  for  each 
case,  we  tried  to  design  a  bonus  system  to  improve  the  quality  of  pyrotechnic  devices. 
We  used  95  %  lower  confidence  bounds  instead  of  the  estimated  overall  reliabilities  to 
decide  whether  to  give  bonuses.  Two  sampling  plans  for  giving  the  bonus  to  the  man¬ 
ufacturers  were  considered: 

■  Single  Sampling  Bonus  System 

■  Double  Sampling  Bonus  System. 

We  simulated  two  sampling  bonus  systems  to  see  the  difference  between  them. 
We  concluded  that  the  double  sampling  system  is  more  sensitive  than  single  sampling 
bonus  system. 

To  formulate  an  approximate  bonus  percentage  as  a  function  of  lower  confidence 
bound  for  bonus  (  LCBFB  ),  the  probability  of  passing  manufacturer  test  Pt,  and  the 
probability  of  passing  the  environment  tests  P2,  we  used  regression  analysis  with 
GRAFSTAT. 
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Because  the  sampling  plan  results  in  so  much  censoring,  the  only  reasonable  esti¬ 
mates  of  overall  reliability  that  we  obtained  were  actually  lower  bounds.  This  makes  the 
bonus  system  conservative  in  the  sense  that  a  bonus  might  not  be  awarded  what  it  is 
deserved.  Therefore,  if  a  bonus  system  is  to  be  implemented,  a  more  comprehensive 
sampling  plan  needs  to  be  devised  which  allows  estimation  of  R.  A  simple  solution  to 
this  problem  can  be  to  apply  all  environmental  tests  to  the  same  sample  which  would 
give  a  measure  of  dependence  between  environmental  tests. 
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APPENDIX  A.  PROGRAM  MLEA 


*  THIS  IS  A  FORTRAN  PROGRAM  TO  CALCULATE  THE  RELIABILITY  OF  AN 

*  ITEM  AFTER  EXPOSURE  TO  SEVERAL  ENVIRONMENTS  WITH  INDEPENDENCE 

*  ASSUMPTION  WHICH  IS  DESCRIBED  IN  CHAPTER  I  .  THE  PROGRAM  READS 

*  5000  SUCCESS  VECTORS,  WHICH  ARE  RANDOMLY  GENERATED  BY  PROGRAM 

*  RANVEC ,  FROM  AN  INPUT  FILE  CALLED  SUCVECT  ONE  BY  ONE.  AFTER 

*  CALCULATION  IT  THEN  PROGRAM  WRITES  ESTIMATED  RELIABILITIES  TO 

*  AN  OUTPUT  FILE  CALLED  RESULT. 


VARIABLES 


k 

SOM 

k 

SOTH 

* 

k 

SOV 

k 

SOA 

k 

R1H 

k 

R2H 

k 

k 

R3H 

k 

R4H 

k 

RHMLE 

k 

k 

N 

k 

X 

k 

FLAG 

k 

R1MAX 

NUMBER  OF  SUCCESSFUL  ITEMS  IN  MANUFACTURER  TEST. 

NUMBER  OF  SUCCESFUL  ITEMS  IN  TEMPERATURE  AND  HUMIDITY 
TEST. 

NUMBER  OF  SUCCESSFUL  ITEMS  IN  VIBRATION  TEST. 

NUMBER  OF  SUCCESSFUL  ITEMS  IN  ALTITUDE  TEST. 

ESTIMATED  PROBABILITY  OF  PASSING  FROM  MANUFACTURER  TEST. 
ESTIMATED  PROBABILITY  OF  PASSING  FROM  TEMPERATURE  AND 
HUMIDITY  TEST. 

ESTIMATED  PROBABILITY  OF  PASSING  FROM  VIBRATION  TEST 
ESTIMATED  PROBABILITY  OF  PASSING  FROM  ALTITUDE  TEST 
ESTIMATED  RELIABILITY  OF  ITEM  AFTER  EXPOSURE  TO 
SEVERAL  ENVIRONMENT  TESTS. 

SAMPLE  SIZES  FOR  EACH  TEST 

NUMBER  OF  SUCCESSFUL  ITEMS  IN  EACH  TEST. 

INDICATOR  VARIABLE  FOR  DETERMINING  EASY  AND  HARD  CASE. 
R1H  VALUE  WHICH  MAXIMIZES  LIKELIHOOD. 
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R2MAX 


* 


* 


* 


* 


R3MAX 

R4MAX 

LMAX 

L 

V 

Q1H 

Q2H 

Q3H 

Q4H 


R2H  VALUE  WHICH  MAXIMIZES  LIKELIHOOD. 

R3H  VALUE  WHICH  MAXIMIZES  LIKELIHOOD. 

R4H  VALUE  WHICH  MAXIMIZES  LIKELIHOOD. 

HIGHEST  MAXIMUM  LIKELIHOOD  VALUE. 

LIKELIHOOD  VALUES  AT  THE  END  OF  EACH  HARD  CASE. 
TOTAL  SAMPLE  SIZES  IN  EACH  HARD  CASE. 

REORGANIZED  PROBABILITY  OF  MANUFACTURER  TEST. 
REORGANIZED  PROBABILITY  OF  TEMP.  AND  HUMIDITY  TEST. 
REORGANIZED  PROBABILITY  OF  VIBRATION  TEST. 
REORGANIZED  PROBABILITY  OF  ALTITUDE  TEST. 


*  TYPE  DECLARATION 

REAL  SOM(5000) ,SOTH(5000) ,SOV(5000) ,SOA(5000) ,R1H ,R2H ,R3H ,R4H , 
+RHMLE(5000) ,N(4) ,X(4) ,R1MAX ,R2MAX ,R3MAX ,R4MAX ,LMAX,V( 7) ,Q1H ,Q2H , 
+Q3H ,Q4H ,NTOT,L( 7 ) 

INTEGER  I,J 
LOGICAL  FLAG(4) 

*  READING  SUCCESS  VECTORS  FROM  SUSVECT  FILE 
DO  60  1=1,5000 

READ( 7 ,*)  SOM(I) ,SOTH  (I) ,SOV(I) ,SOA(I) 

*  rt.'.-A-  A  *A-A  it  *HWWHHWHHWWHW******  A  A' A  A  A  is  ft***  it  A-A-A  A  A  .'.-AA  A-A-.'-AA 

*  FILES  FOR  READING  AND  WRITING 

CALL  EXCMS  ('FILEDEF  7  DISK  SUCVECT  DATA  Al') 

CALL  EXCMS  (’FILEDEF  16  DISK  RESULT  DATA  Al’) 

A  ** A  A  it**i;it*it**ftirieic*-h~/t-ie**-icA"is  A  A  A  A  A  irMrMrMt*******  Me**  A  A*  A-  A  A  **  A-  *****  A 

*  INITIALIZATION  OF  SAMPLE  SIZES 

N( 1 )=  20.0 

N ( 2 ) =  20.0 
N(3)=  32.0 
N(4)=  20.0 

*  INITIALIZATION  OF  FLAG  VARIABLES 

DO  10  J=2 ,4 

FLAG( J)=  .FALSE. 
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10  CONTINUE 


NUMBER  OF  SUCCESSES  IN  EACH  TEST 
X( 1 )=  SOM(I) 

X(2)=  SOTH(I) 

X( 3)=  SOV(I) 

X(4)=  SOA(I) 


CHECK  OPERATION  FOR  EASY  CASE 
IF  ( (X( l)/20 . 00) .GE.(X(2)/20.00))  THEN 
FLAG( 2)=  .TRUE. 

END  IF 

IF  ((X(l)/20.00).GE.(X(3)/32.00))  THEN 
FLAG(3)=  .TRUE. 

END  IF 

IF  ( (X( 1 )/20 . 00) ,GE . (X(4)/20 . 00) )  THEN 
FLAG(4)=  .TRUE. 

END  IF 

IF  (  FLAG(2).AND.FLAG(3).AND.FLAG(4)  )  THEN 


CALCULATIONS  IN  EASY  CASE 
R1H=  X(l)/20.0 
R2H=  X(2)/20 . 0 
R3H=  X(3)/32.0 
R4H=  X(4)/20 . 0 
Q1H=  R1H 
Q2H=  R2H/Q1H 
Q3H=  R3H/Q1H 
Q4H=  R4H/Q1H 

RHMLE( I )=  Q1H*Q2H*Q3H*Q4H 
GO  TO  50 
END  IF 


CALCULATIONS  IN  HARD  CASES 
INITIALIZATION 
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DO  30  J=l,7 
V(J)=  0.0 
L( J)=  0.0 
30  CONTINUE 
LMAX=  0.0 
NTOT=  0.0 


40 


CASE  1  R1H.LT.1.0  AND  R1H  =  R2H  =  R3H  =  R4H 

DO  40  J=1 , 4 

V(l)=  V(1)+X(J) 

NTOT=  NTOT+N(J) 

CONTINUE 
R1H=  V(l)/NTOT 

L( 1 )=  ( (R1H)**V( 1) )*( (NTOT-V( 1) ) /NTOT)**(NTOT-V( 1) ) 
IF  (L(l).GT.LMAX)  THEN 
LMAX=  L( 1 ) 

R1MAX=  R1H 
R2MAX=  R1H 
R3MAX=  R1H 
R4MAX=  R1H 
END  IF 


CASE  2  RIH. LT. 1.0,  RIH  =  R2H  =  R3H  ,  R4H  IS  BETWEEN  0.0  AND  R1H 
IF  (  X(4)/20 . 0 . NE . 1 . 0  )  THEN 
V(2)=  X(l)+X(2)+X(3) 

NTOT=  N(I)+N(2)+N(3) 

R1H=  V(2)/NTOT 
IF((X(4)/20.00).LE.R1H)  THEN 

L(2)=  ( (R1H)**V(2) )*( (NTOT-V( 2) )/NTOT)**(NTOT-V(2) )*( (X(4) 
+/20.0)**X(4))*(((20.0-X(4))/20.0)**(20.0-X(4))) 

IF  (L(2) .GT.LMAX)  THEN 
LMAX=  L( 2) 

R1MAX=  RIH 
R2MAX=  RIH 
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R3MAX=  R1H 
R4MAX=  X(4)/20.0 
END  IF 
END  IF 
END  IF 


*  CASE  3  R1H.LT.1.0,  R1H  =  R2H  =  R4H  ,  R3H  IS  BETWEEN  0.0  AND  R1H 
IF  (  X( 3)/32 . 0 . NE . 1 . 0  )  THEN 

V(3)=  X(l)+X(2)+X(4) 

NTOT=  N(l)+N(2)+N(4) 

R1H=  V(3)/NTOT 

IF  ((X(3)/32.00) .LE.R1H)  THEN 

L(3)=  ( (R1H)**V( 3) )*( (NTOT~V( 3) )/NTOT)**(NTOT-V(3) )*( (X( 3) 
+/32 . 0)**X( 3) )*( ( (32 . 0-X(3) )/32 . 0)**(32 . 0-X(3) ) ) 

IF  (L(3) .GT.LMAX)  THEN 
LMAX=  L(3) 

R1MAX=  R1H 
R2MAX=  R1H 
R3MAX=  X( 3)/32 . 0 
R4MAX=  R1H 
END  IF 
END  IF 
END  IF 

*  CASE  4  R1H.LT. 1.0,  R1H  =  R3H  =  R4H  ,  R2H  IS  BETWEEN  0.0  AND  R1H 
IF  (  X(2)/20.0.NE. 1.0  )  THEN 

V(4)=  X(l)+X(3)+X(4) 

NTOT=  N(l)+N(3)+N(4) 

R1H=  V(4)/NT0T 

IF  ((X(2)/20.00) .LE.R1H)  THEN 

L(4)=  ( (R1H)**V(4) )*( (NT0T-V(4) )/NT0T)**(NT0T-V(4) )*((X(2) 
+/20.0)**X(2))*(((20.0-X(2))/20.0)**(20.0-X(2))) 

IF  (L(4). GT.LMAX)  THEN 
LMAX=  L(4) 
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* 


R1MAX=  R1H 
R2MAX=  X(2)/20 . 0 
R3MAX=  R1H 
R4MAX=  R1H 
END  IF 
END  IF 
END  IF 


CASE  5  R1H.LT.1.0,  R1H  =  R2H  R3H  ,  R4H  ARE  BETWEEN  0.0  AND  R1H 
IF  ((X(3)/32.0.NE. 1.0) .AND. (X(4)/20 . 0 . NE . 1.0))  THEN 
V(5)=  X(l)+X(2) 

NTOT=  N(l)+N(2) 

R1H=  V(5)/NTOT 

IF((X(3)/32.00).LE . R1H . AND . (X(4)/20 . 00) .LE . R1H)  THEN 

L(5)=  ( (R1H)**V( 5) )*( (NTOT'V(5) )/NTOT)**(NTOT-V(5 ) )*((X(3) 
+/32.0)**X(3))*((32.0-X(3))/32.0)**(32.0-X(3))*((X(4)/20.0)** 
+x(4) )*( (20 . 0-X(4) )/20 .0)**(20 . 0-X(4) ) 

IF  (L(5) .GT.LMAX)  THEN 
LMAX=  L(5) 

R1MAX=  R1H 
R2MAX=  R1H 
R3MAX=  X( 3)/32 . 0 
R4MAX=  X(4)/20 . 0 
END  IF 
END  IF 
END  IF 


*  CASE  6  R1H.LT. 1 .0 ,  R1H  =  R3H  R2H  ,  R4H  ARE  BETWEEN  0.0  AND  R1H 
IF  ((X(2)/20.0.NE.1.0) .AND ,(X(4)/20.0.NE.1.0))  THEN 
V(6)=  X(l)+X(3) 

NTOT=  N( 1)+N( 3) 

R1H=  V(6)/NTOT 

IF((X(2)/20.00). LE ,R1H . AND . (X(4)/20 . 00) . LE .R1H)  THEN 
L(6)=  ( (R1H)**V(6) )*( (NTOT-V(6) )/NTOT)**(NTOT-V(6) )*( (X(2) 
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+/20.0)**X(2))*((20.0-X(2))/20.0)**(20.0-X(2))*((X(4)/20.0)** 
+x(4) )*( (20 . 0-X(4) )/20 . 0)**(20 .0-X(4) ) 

IF  (L(6) .GT.LMAX)  THEN 
LMAX=  L(6) 

R 1MAX=  R1H 
R2MAX=  X(2)/20 . 0 
R3MAX=  R1H 
R4MAX=  X(4)/20.0 
END  IF 
END  IF 
END  IF 


CASE  7  R1H.LT.1.0,  R1H  =  R4H  R2H  ,  R3H  ARE  BETWEEN  0.0  AND  R1H 
IF  ((X(2)/20.0.NE. 1.0) . AND. (X(3)/32.0.NE.1.0))  THEN 
V( 7 ) =  X(l)+X(4) 

NTOT=  N(l)+N(4) 

R1H=  V(7)/NTOT 

IF( (X(2)/20. 00). LE.R1H.AND.(X( 3)/32. 00). LE.R1H)  THEN 
L(7)=  ((RlH)**V(7))*((NTOT-V(7))/NTOT)**(NTOT-V(7))*((X(3) 
+/32.0)**X(3))*((32.0-X(3))/32.0)**(32.0-X(3))*((X(2)/20.0)**X(2)) 
+*((20.0-X(2))/20.0)**(20.0-X(2)) 

IF  (L(7) .GT.LMAX)  THEN 
LMAX=  L(7) 

R1MAX=  R1H 
R2MAX=  X(2)/20 .0 
R3MAX=  X(3)/32.0 
R4MAX=  R1H 


CALCULATION  OF  RHMLE  VALUE  WITH  RESPECT  TO  CASE  WHICH  HAS  LARGEST 

MAXIMUM  LIKELIHOOD  IN  HARD  CASE 

R1H=  R1MAX 

R2H=  R2MAX 

R3H=  R3MAX 

R4H=  R4MAX 

Q1H=  R1H 

Q2H=  R2H/Q1H 

Q3H=  R3H/Q1H 

Q4H=  R4H/Q1H 

RHMLE( I )=  Q1H*Q2H*Q3H*Q4H 


****** 


^  ^  »  ‘,V  » •fv*  "A*  »•  Vr  A  tV  » ‘A‘  V  ".V  t1:  »  »  Vr  »  »  Vr 


WRITING  AFTER  EACH  CALCULATION 
50  WRITE  (16,*)  RHMLE ( I ) 

60  CONTINUE 
STOP 
END 
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APPENDIX  B.  PROGRAM  RANVEC 


PROGRAM  RANVEC 


* 


THIS  IS  THE  PROGRAM  TO  GENARETE  RANDOM  NUMBERS  FROM  BINOMIAL 
DISTRIBUTION. THE  PROGRAM  READS  PROBABILITIES  OF  BEING  SUCCESSFUL 
IN  FOUR  TESTS  INTERACTIVELY.  IT  GENERATES  UNIFORMLY  DISTRIBUTED 
RANDOM  NUMBERS  WITH  THESE  PROBABILITIES  ACORDING  TO  SAMPLE  SIZE 
OF  EACH  TEST.  FOR  EACH  TEST  PROGRAM  COUNTS  UNIFORMLY  DISTRIBUTED 
RANDOM  NUMBERS,  WHICH  HAVE  GREATER  THAN  OR  EQUAL  PROBABILITY 
WITH  RESPECT  TO  THE  GIVEN  PROBABILITY  FOR  THAT  TEST. TOTAL  COUNTS 
GIVE  US  SUCCESSFUL  ITEM  NUMBERS  FOR  EACH  TEST.  PROGRAM  UPDATES 
SEEDS  AND  CALLS  SUBROUTINE  RANNUM  IN  EACH  ITERATION.  THE  PROGRAM 
GENERATES  5000  SUCCESS  VECTOR  AND  WRITES  THEM  TO  AN  OUTPUT  FILE 
CALLED  SUSVECT. 


kkkkkkkkkkk^rf.-k-frkkkkkk-i 


VARIABLES 


k 

k 

k 


* 


*• 


k 


k 


k 


k 


k 


k 


PSIM  :  PROBABILITY  OF  SUCCESS  IN  MANUFACTURER  TEST. 

PSITH  :  PROBABILITY  OF  SUCCESS  IN  TEMP.  AND  HUMIDITY  TEST. 

PSIV  :  PROBABILITY  OF  SUCCESS  IN  VIBRATION  TEST 

PSIA  :  PROBABILITY  OF  SUCCESS  IN  ALTITUDE  TEST 

NUM1  :  COUNTER  FOR  MANUFACTURER  TEST 

NUM2  :  COUNTER  FOR  TEMP.  AND  HUMIDITY  TEST 

NUM3  :  COUNTER  FOR  VIBRATION  TEST 

NUM4  :  COUNTER  FOR  ALTIUTDE  TEST 

X  :  BINOMIALLY  DISTRIBUTED  RANDOM  NUMBER  FOR  MANUFACTURER 

TEST 

V  :  BINOMIALLY  DISTRIBUTED  RANDOM  NUMBER  FOR  TEMP. AND 

HUMUDITY  TEST 
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Y 

Z 

A 

B 

C 

D 

ISEED 

KSEED 

LSEED 

MSEED 


BINOMIALLY  DISTRIBUTED  RANDOM  NUMBER  FOR  VIBRATION  TEST 
BINOMIALLY  DISTRIBUTED  RANDOM  NUMBER  FOR  ALTITUDE  TEST 
UNIFORMLY  DISTRIBUTED  RANDOM  NUMBER  FOR  MANUFACTURER 
TEST 

UNIFORMLY  DISTRIBUTED  RANDOM  NUMBER  FOR  TEMP. AND 
HUMIDITY  TEST 

UNIFORMLY  DISTRIBUTED  RANDOM  NUMBER  FOR  VIBRATION  TEST 

UNIFORMLY  DISTRIBUTED  RANDOM  NUMBER  FOR  ALTITUDE  TEST 

SEED  NUMBER  FOR  MANUFACTURER  TEST 

SEED  NUMBER  FOR  TEMP.  AND  HUMIDITY  TEST 

SEED  NUMBER  FOR  VIBRATION  TEST 

SEED  NUMBER  FOR  ALTITUDE  TEST 


t*****-i 


«  A  <V  A  A A  A  A  A  /, 


******** 


TYPE  DECLARATION 

REAL  PSIM , PSITH ,PSIV, PSIA ,NUM1 ,NUM2 ,NUM3 ,NUM4 ,X(5000) ,  V(5000) 
+,Y(5000) ,Z(5000) ,A,B,C,D 
INTEGER  I SEED, KSEED, LSEED, MSEED 

ftftftft'.Vftft  /■■■■’.  *  A  A  ft  A  A  A*A-ftftftA****ft***ft****ft*ft**ft*-ft-*ft**ft*****ftft  A  A  A  ***** 

INITIALIZATION 
ISEED  =  45267 
KSEED  =  113234 
LSEED  =  435 
MSEED  =  1 

ft***********  A  AAA  **************  **********  ***  **********  A  ft  A  A  A  A  ***** 

READING  TEST  PROBABILITIES  OF  BEING  SUCCESSFUL  IN  EACH  TEST 
WRITEC*,*)' PLEASE  WRITE  PROBABILITY  OF  BEING  SUCCESSFUL  IN 
+  MANUFACTURER  TEST' 

READ  (*,*)  PSIM 

WRITE(*,*) 'PLEASE  WRITE  PROBABILITY  OF  BEING  SUCCESSFUL  IN 
+  TEMPERATURE  AND  HUMIDITY  TEST' 

READ  (*,*)  PSITH 

WRITEC*,*) 'PLEASE  WRITE  PROBABILITY  OF  BEING  SUCCESSFUL  IN 
+  VIBRATION  TEST’ 

READ  C*,*)  PSIV 


WRITE(*,*) 'PLEASE  WRITE  PROBABILITY  OF  BEING  SUCCESSFUL  IN 
+  ALTITIDE  TEST' 

READ  (*,*)  PSIA 


JLJUJUL^JLJL^LJUAJLJUUJUJULJULJbJLj^kJL^LJLJLJLJUJUJU^JLJ^JLJLJUJLJULJUJUJU^^^U^UUJ^^UJ^LJUAJ^LJUJkJLJUUX 

OPENING  AN  OUTPUT  FILE  TO  WRITE  THE  RESULTS 
CALL  EXCMS  ( 'FILEDEF  16  DISK  SUSVECT  DATA  Al') 


GENERATION 
DO  50  J=l, 5000,1 


INITIALIZATION  IN  EACH  ITERATION 
NUM1=  0.0 


NUM2=  0.0 
NUM3=  0.0 


NUM4=  0.0 


SEEDS  UPDATATION  IN  EACH  ITERATION 
ISEED=ISEED+17 
KSEED=KSEED+1356 
LSEED=LSEED+1 
MSEED=MSEED+789 


**  /.  A  A  A'  A1  ft  /. A-*-.1.-*  A  ***********  ft  ft  A  ft  A-ft rft-ft-ft-A-ft  ft  ft  ft  *********  A  ftft-ft-ftftrft-ftftrft-ftftT'cftr-ft 

GENERATION  OF  SUCCESSFUL  NUMBER  OF  ITEMS  IN  MANUFACTURER  TEST 
DO  10  1=1,20,1 

CALL  RANNUM  (1,1 SEED ,0.0, 1.0, 0.0, A) 

IF  (  A.LT.PSIM  )  THEN 
NUM1=NUM1+1 
END  IF 
X(J)  =  NUM1 
10  CONTINUE 


GENERATION  OF  SUCCESSFUL  NUMBER  OF  ITEMS  IN  TEMP.  AND  HUM.  TEST 
DO  20  K=1 ,20,1 

CALL  RANNUM  ( 1 ,KSEED ,0 . 0 , 1 . 0 ,0 . 0 ,B) 
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20 


IF  (  B.LT.PSITH  )  THEN 
NUM2=NUM2+1 
END  IF 
V(J)  =  NUM2 
CONTINUE 


30 


GENERATION  OF  SUCCESSFUL  NUMBER  OF  ITEMS  IN  VIBRATION  TEST 
DO  30  L=1 , 32 , 1 

CALL  RANNUM  ( 1  ,LSEED,0 . 0 , 1 . 0 , 0 . 0 ,C) 

IF  (  C.LT.PSIV  )  THEN 
NUM3=NUM3+1 
END  IF 
Y(J)  =  NUM3 
CONTINUE 


40 


GENERATION  OF  SUCCESSFUL  NUMBER  OF  ITEMS  IN  ALTITUDE  TEST 
DO  40  M=1 ,20,1 

CALL  RANNUM  ( 1  ,MSEED,0 . 0 , 1 . 0 ,0 . 0 ,D) 

IF  (  D.LT.PSIA  )  THEN 
NUM4=NUM4+1 
END  IF 
Z(J)  =  NUM4 
CONTINUE 


WRITING  THE  RESULTS  TO  AN  OUTPUT  FILE  AS  4  TUPLE 
WRITE  (16,1)  X(J),V(J),Y(J),Z(J) 

1  FORMAT  ( IX ,F12 . 7 ,4X ,F12 . 7 ,4X,F12 . 7 ,4X,F12 . 7 ,4X) 

50  CONTINUE 
STOP 
END 

/■  A".V'A"A’-A-A-A-A  A  A  A  A  A  'h  iV  A1  !,~K  A  it  Vi  ]V  i;  -,'r  -k  i?  it f-k-K-jrh'  A  A 

SUBROUTINE  RANNUM (DISTN ,  SEED,  RPARM1 ,  RPARM2 ,  IPARM,  X) 

•frtrMrfrfc frfc/ck  A’  it  it  ‘A*  ic  fr  ic  i:  iciv  ir/rfn A’  it  IcJrtrfc  is  is  -frfrfrfrfrfHr Jcirirts  is^rk^sitirkr/riHc^ririsis'is'isisisisis'isis 

THIS  SUBROUTINE  IS  A  PART  OF  SIMUTIL  FORTRAN  WHICH  IS  WRITTEN 
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BY  DR.  M.  P.  BAILEY.  THIS  SUBROUTINE  PROVIDES  AN  INTERFACE  WITH 
THE  LLRANDOMII  ROUTINES  PROVIDED  IN  THE  NONIMSL  LIBRARY.  THE 
PARAMETER  REQUIRMENTS  AND  CALLING  PROCEDURES  ARE  AS  FOLLOWS: 
DISTN  =  DISTRIBUTION  TYPE  YOU  WANT  TO  SELECT  AN  INTEGER  BETWEEN  1 
AND  7. 

SEED  =  THE  RANDOM  NUMBER  SEED  YOU  WISH  TO  USE. 

RPARM1 ,  RPARM2 ,  AND  IPARM  ARE  REAL  AND  INTEGER  PARAMETERS  PASSED 
TO  THE  ROUTINE  WITH  MEANINGS  WHICH  VARY  WITH  THE  TYPE  OF  DISTRI_ 
BUTION  YOU  DESIRE. 

X  =  THE  RETURNED  RANDOM  NUMBER,  IT  IS  ALWAYS  REAL. 

DISTRIBUTION  NUMBERS  AND  THE  ASSOCIATED  PARM  DEFINITIONS 

1-  -UNIFORM  ON  THE  INTERVAL  RPARM1  TO  R PARM 2 . 

2- -NORMAL  WITH  MEAN  RPARM1  AND  VARIANCE  RPARM2 . 

3- -EXPONENTIAL  WITH  RATE  RPARM1 . 

4- -C0UCHY  WITH  A  =  RPARM1  AND  B  =  RPARM2 . 

5-  -GAMMA  WITH  SHAPE  R PARM 2  AND  RATE  RPARM1 . 

6-  -POISSON  WITH  RATE  RPARM1 . 

7-  -GEOMETRIC  WITH  P  =  RPARM1 . 

k-kkkirk'frkkkk-k-kkirk-krkkkrkkk-k-k-k  kk-kkk  kkkkk-kk-kkkk  k  kk  k-kkk-kk-.‘.~k-k-kk-k-k-k-kk-k-k-k 

TYPE  DECLARATION 

REAL  RPARM 1 , RPARM2 , X , TEMP , VARI AT( 1 ) 

INTEGER  DISTN,  SEED,  IPARM,  N 

IF  (DISTN. LE.O. OR. DISTN. GT. 8)  THEN 

WRITEUO,  *)  'ILLEGAL  CALL  TO  RANNUM,  BAD  DISTN' 

STOP 

ENDIF 

GOTO  (10,  20,  30,  40,  50,  60,  70),  DISTN 

kkk A  k  k  k  k  k  A  /.  k  /.  kk  kk-k-A  k * k  /<■  kkie-k  A  A  A  A  A  irkirkkkkkkkkkkkkkkk-kkkkkk-kkk-kk-kk 


10 


GENERATE  A  UNIFORM  BETWEEN  RPARM1  AND  RPARM2 
CONTINUE 

IF  (RPARM I  -  RPARM2.EQ.0)  THEN 

WRITEUO,  *)  'ILLEGAL  EQUAL  RPARMS  IN  RANNUM' 
STOP 
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ENDIF 

IF  ( RPARM 1 . GT . RPARM2 )  THEN 
TEMP  =  RPARM 1 
RPARM 1  =  RPARM2 
RPARM2  =  TEMP 
ENDIF 

CALL  LRND(SEED,  VARIAT,  1,  1,  0) 

VARIAT(l)  =  RPARM 1  +  (RPARM2  -  RPARM1)  *  VARIAT(l) 
GOTO  80 


GENERATE  A  NORMAL  WITH  MEAN  RPARM1  AND  STDDEV  RPARM2 
CALL  LNORM(SEED,  VARIAT,  1,  1,  0) 

VARIAT(l)  =  (VARIAT(l)  *  RPARM2)  +  RPARM 1 
GOTO  80 


JUJU 


GENERATE  AN  EXPONENTIAL  WITH  RATE  (1/MEAN)  RPARM 1 
CONTINUE 

IF  (RPARM 1.EQ.0)  THEN 

WRITE(10,  *)  'ILLEGAL  ZERO  RATE  IN  RANNUM' 

STOP 

ENDIF 

CALL  LEXPN(SEED,  VARIAT,  1,  1,  0) 

VARIAT(l)  =  VARIAT(l)  /  RPARM 1 
GOTO  80 

A- A-  i;  A-  A  -A -A-  A  A-  A  A-AA-AArA-A  AAA  AAAAA  'A  A  'A'  A  'A"A"A  A  A'A'A'A'A'A'A  A  A  A  A  A  A  A  AA  A  A  A- A~A-A  A- A- A  A-  .cAAA 

GENERATE  A  COUCHY  WITH  A  =  RPARM 1  AND  B  =  RPARM 2 
CONTINUE 

IF  (RPARM2.LE.0)  THEN 

WRITE(10,  *)  'ILLEGAL  COUCHY  SPREAD  IN  RANNUM,  B  =  ’ ,RPARM2 
STOP 
ENDIF 

CALL  LCCHY(SEED,  VARIAT,  1,  1,  0) 

VARIAT(l)  =  (VARIAT(l)  *  RPARM2)  +  RPARM1 


GOTO  80 


*  GENERATE  GAMMA  WITH  SHAPE  RPARM'  AND  RATE  RPRAM1 
50  CONTINUE 

IF  ( RPARM l.LE.O)  THEN 

WRITE(10,  *)  'ILLEGAL  NONPOSITIVE  GAMMA  RATE  IN  RANNUM’ 
STOP 
END  IF 

IF  (RPARM2.LE.0)  THEN 

WRITE( 10 ,  *)  'ILLEGAL  SHAPE  PARAMETER  IN  RANNUM' 

STOP 

ENDIF 

CALL  LGAMA(SEED,  VARIAT,  1,  1,  0,  RPARM2) 

VARIAT(l)  =  VARIAT(l)  *  (1.0  /  RPARM 1 ) 

GOTO  80 

*  GENERATE  POISSON  WITH  RATE  RPRAM1 
60  CONTINUE 

IF  (RPARM l.LE.O)  THEN 

WRITE(10,  *)  'ILLEGAL  POISSON  RATE  IN  RANNUM' 

STOP 

ENDIF 

CALL  LPOIS(SEED,  VARIAT,  1,  1,  0,  RPARM 1) 

GOTO  80 

*  GENERATE  GEOMETRIC  WITH  P  =  RPRAM1 
70  CONTINUE 

IF  (RPARM l.LE.O)  THEN 

WRITE(10,  *)  'ILLEGAL  GEOM  PROB  IN  RANNUM' 


ENDIF 

CALL  LGEOM(SEED,  VARIAT,  1,  1,  0,  RPARM 1) 
GOTO  80 


80  CONTINUE 


91 


X  =  VARIAT(l) 
END 
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APPENDIX  C.  PROGRAM  SORT 


THIS  IS  THE  SORTING  PROGRAM.  PROGRAM  USES  BUBBLE  SORT  ALGORITHM. 
PROGRAM  READS  ESTIMATED  RELIABILITIES  FROM  AN  INPUT  FILE  CALLED 
RESULT.  IT  SORTS  FROM  SMALLEST  TO  LARGEST,  AND  WRITES  IN  TO  AN 
OUTPUT  FILE  CALLED  FRESULT  WITH  95  %  LOWER  CONFIDENCE  BOUND. 


* 


VARIABLES 


A  :  ESTIMATED  RELIABILITY 

FLAG  :  INDICATOR  VARIABLE  TAKES  VALUE  '  OF  '  AND  '  OFF  ' 

^XXAXA.»t.J,JL;)lt.^.XJU^.L^L^J.^LjL^^.I.AJL^LjLjL^LX^JLJLJLXJLJL.X^L;L^LJL..JL^JkXX:LJLXA^XJUJLjUA;LXJLXJUJL 


TYPE  DECLARATION 
CHARACTER  FLAG* 3 
REAL  A(5000) 

INTEGER  I ,N, J 

OPENING  AN  INPUT  AND  AN  OUTPUT  FILE 

CALL  EXCMS  ('FILZDEF  9  DISK  RESULT  DATA  Al') 

CALL  EXCMS  ( ' FILEDEF  15  DISK  FRESULT  DATA  Al') 

A  h  A’  A  A  A' A  A  i\ 'A' »  vV tV Vr* VrVr •k-k'k-k-lr'fc'k'tr-k'k'trix VnVVr -fr 


READING  ESTIMATED  RELIABILITIES 
DO  10  1=1,5000 
READ(9 ,*)  A( I ) 

10  CONTINUE 


SORTING  OPERATION 
N= I  - 1 
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DO  30  I=N ,2,-1 
FLAG= ' OFF 1 


DO  20  J=1,I-1 

IF  (A(J).GT.A(J+1))  THEN 
TEMP=A( J) 

A( J)=A( J+l) 
A(J+1)=TEMP 
FLAG= ' ON ' 

END  IF 

20  CONTINUE 

IF  ( FLAG. EQ. 'OFF')  THEN 
GO  TO  40 
END  IF 

30  CONTINUE 
40  CONTINUE 


WRITING  THE  RESULTS  IN  ASCENDING  ORDER  TO  OUTPUT  FILE 
DO  50  1=1,5000 

WRITE  (15,*)  A(I) 

50  CONTINUE 

WRITE  (15,1)  A(250) 

1  FORMAT  (/// , 15X, ' 95  J,  LOWER  CONFIDENCE  BOUND  IS’ , 1X,F12. 7) 
STOP 
END 
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APPENDIX  D.  PROGRAM  INITIAL 


PROGRAM  INITIAL 


THIS  PROGRAM,  CALCULATES  INITIAL  GUESSES  FOR  PARAMETERS  IN 
LOGLINEAR  MODEL  BY  MEANS  OF  PROGRAM  PARAM  WHICH  IS  IN  APPENDIX  E. 
PROGRAM  SUPPLIES  PARTIAL  SUMS  OF  EXPECTATION,  TO  SOLVE  EQUATIONS 
IN  PROGRAM  PARAM.  IT  READS  INTERACTIVELY  NUMBER  OF  FAILURES  IN 
TESTS.  PROGRAM  WRITES  RESULTS  TO  AN  OUTPUT  FILE  CALLED  EXPECT. 


VARIABLES 


FOM 

FOTH 

FOV 

FOA 

PI 

P2 

P3 

P4 

Q1 

Q2 

Q3 

Q4 

XI 

X2 

X3 

X4 


NUMBER  OF  FAILURES  IN  MANUFACTURER  TEST 

NUMBER  OF  FAILURES  IN  TEMPERATURE  AND  HUMIDITY  TEST 

NUMBER  OF  FAILURES  IN  VIBRATION  TEST 

NUMBER  OF  FAILURES  IN  ALTITUDE  TEST 

SUCCESS  RATIO  FOR  MANUFACTURER  TEST 

SUCCESS  RATIO  FOR  TEMPERATURE  AND  HUMIDITY  TEST 

SUCCESS  RATIO  FOR  VIBRATION  TEST 

SUCCESS  RATIO  FOR  ALTITUDE  TEST 

FAILURE  RATIO  FOR  MANUFACTURER  TEST 

FAILURE  RATIO  FOR  TEMPERATURE  AND  HUMIDITY  TEST 

FAILURE  RATIO  FOR  VIBRATION  TEST 

FAILURE  RATIO  FOR  ALTITUDE  TEST 

PARTIAL  SUM  OF  EXPECTED  NUMBERS  IN  CELLS  FOR  EQUATION  1 

PARTIAL  SUM  OF  EXPECTED  NUMBERS  IN  CELLS  FOR  EQUATION  2 

PARTIAL  SUM  OF  EXPECTED  NUMBERS  IN  CELLS  FOR  EQUATION  3 

PARTIAL  SUM  OF  EXPECTED  NUMBERS  IN  CELLS  FOR  EQUATION  4 


TYPE  DECLARATION 
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REAL  FOM ,FOTH,FOV ,FOA ,P1 ,P2 ,P3 ,P4 ,Q1 ,Q2 ,Q3 ,Q4 ,PROD,Xl ,X2 ,X3 ,X4 


*  OPENING  AN  OUTPUT  FILE 

CALL  EXCMS  ( 'FILEDEF  13  DISK  EXPECT  DATA  Al') 


READINC  NUMBER  OF  FAILURES  IN  EACH  TEST  INTERACTIVELY 
WRITE!*,*) 'PLEASE  WRITE  NUMBER  OF  FAILURES  IN  MANUFACTURER  TEST' 
READ(* ,*)  FOM 

WRITE!*,*) 'PLEASE  WRITE  NUMBER  OF  FAILURES  IN  TEMP. AND  HUM. TEST’ 
READ!* ,*)  FOTH 

WRITE!*,*) 'PLEASE  WRITE  NUMBER  OF  FAILURES  IN  VIBRATION  TEST' 
READ!*,*)  FOV 

WRITE(*,*) 'PLEASE  WRITE  NUMBER  OF  FAILURES  IN  ALTITUDE  TEST' 
READ(*,*)  FOA 


CALACULATION  OF  SUCCESS  RATIOS 

Pl=( 20 . O-FOM) /8 . 00 

P2=(20 . 0-FOTH)/8 . 00 

P3=( 32 . 0-F0V)/8 . 00 

P4=( 20 . 0-F0A)/8 . 00 

CALCULATION  OF  FAILURE  RATIOS 

Ql=F0M/8 . 00 

Q2=F0TH/8 . 00 

Q3=F0V/8 . 00 

Q4=F0A/8 . 00 


*  CALCULATION  PARTIAL  SUMS  OF  EXPECTATION  FOR  PROGRAM  PARAM 
PR0D=2 . 0*P3+2 . 0*Q3+2 . 0*P4+2 . 0*Q4 
X 1=4 . 0*P 1+4 . 0*P2+PROD 
X2=4 . 0*Pl+4 . 0*Q2+PROD 
X3=4 . 0*Ql+4 . 0*P2+PROD 
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X4=4 . 0*Q 1+4 . 0*Q2+PROD 


WRITING  THE  RESULTS  TO  AN  OUTPUT  FILE 

WRITE  (13,*)  X1.X2.X3.X4 

STOP 

END 
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APPENDIX  E.  PROGRAM  PARAM 


PROGRAM  PARAM 


THIS  IS  THE  PROGRAM  TO  CALCULATE  PARAMETERS  OF  LOGLINEAR  MODEL 
WITH  IMSL  SUBROUTINE.  IT  TAKES  PARTIAL  SUMS  FROM  PROGRAM  INITIAL 
AND  SOLVES  FOUR  NONLINEAR  EQUATIONS,  WHICH  HAVE  FOUR  UNKNOWNS. 
THE  PROGRAM  USES  AN  IMSL  SUBROUTINE  CALLED  DNEQNF  TO  SOLVE  THIS 
EQUATION.IT  WRITES  SOLUTIONS  OF  EQUATIONS  TO  OUTPUT  FILE  CALLED 
PARAM  DATA. 


VARIABLES 


ITMAX  :  MAXIMUM  ITERATION  NUMBER. 

N  :  PARAMETER 

XGUESS  :  INITIAL  GUESS  FOR  FOUR  NONLINEAR  EQUATIONS. 

F  :  NONLINEAR  EQUATIONS. 

•frfrfrfcirfrfrfHHc A  />  A  &  A  -kk  A  k  fr/cJcfcjcie'to  A1  ‘A’  A  A  A-A-tV  irk  •/: ******* VrA”A-  A  rfn  *  A"/ "A1  *■*■■*■  VrA- 

TYPE  DECLARATION 


PARAMETER  (N=4) 

REAL*8  ERRREL 

* 

INTEGER  ITMAX, N 

M 

INTEGER  K 

REAL*8  FNORM,X(N),XGUESS(N) 

EXTERNAL  ACN 

jl.  j.j  jl  j  j,  j. ■  i  t  J,  ,r  •*.-  A-A- -A- -A-.W-  -A-  -A  -A- A-  A-  A  -A--A- -A-  -A-Ar  -A-A 

*  OPENING  A  FILE  FOR  WRITING  RESULTS 
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CALL  EXCMS  ( 'FILEDEF  9  DISK  PARAM  DATA  A') 


INITIAL  GUESS 

DATA  XGUESS/3 .500,3.500,3.500,3. 5D0/ 


INITIALIZATION 
ERRREL  =  0 . 000  IDO 
ITMAX  =  10000 


CALLING  OF  IMSL  SUBROUTINE 

CALL  DNEQNF  (  ACN, ERRREL, N, ITMAX, XGUESS ,X ,FNORM) 


RESULTS 

WRITE  (9,1)  (X(K),K=l,N),FNORM 

FORMAT('THE  SOLUTION  TO  THE  SYSTEM  IS X=( ' ,4F8 . 2 WITH 
+FNORM= ' , F8 . 2 , // ) 

END 


SUBROUTINE  ACN  (X,F,N) 

kkkkkkkkkkkkkkkkkk-M.-kkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkk 

VARIABLES 

X  :  INITIAL  GUESS 
F  :  NONLINEAR  EQUATIONS 

kk  kkkkkkkkkkkkkk  kkkiik-frfckkkkkkkkkkkkkkkkkkk  kkkkkkkkkkkkkkkkkkkkkk 


TYPE  DECLARATION 
REAL*8  X(N) ,F(N) 
INTEGER  N 


«V 'A* ’/?»* 1 A* *« >* 'A* '<< ‘A* * A* *A*  A  A  /%  A  A1  'A“<V  mJ%  A1  »  A* 'A* *A* «»* *»V /> ^ «V 'A1 » '/>  <V *« ** '/>  A **c A «V«»'«VV>'«\,Vo«ViV*A*V\'A'vV'/'>*A'^ 


1  ST  EQUATION 

F(1)=DEXP(X(1))*(DEXP(X(2)+X(3)+X(4))+ 

+  (1/DEXP(X(2)+X(3)+X(4)))+ 

+  (DEXP(X(2)+X(3) )/DEXP(X(4) ) )+ 

+  (DEXP(X(2)+X(4) )/DEXP(X( 3) ) )+(DEXP(X( 3) 

+  +X(4))/DEXP(X(2)))+(DEXP(X(2))/DEXP(X(3)+X(4)))+ 


+  (DEXP(X(3))/DEXP(X(2) 

+  +X( 3) ) )+(DEXP(X(4) )/DEXP(X(2)+X( 3) ) ) ) -LOG (32 .00) 


t************* 


r**********************- 


t***** 


2  ST  EQUATION 

F( 2)=DEXP(X( 1 )+X(2) )*(DEXP(X( 3)+X(4) )+ 

+  ( 1/DEXP(X( 3)+X(4) ) )+(DEXP(X( 3) )/ 

+  DEXP(X(4) ) )+(DEXP(X(4) ) /DEXP(X( 3) ) ) ) -LOG (23 .00) 


3  ST  EQUATION 

F( 3)=DEXP(X( 1)+X( 3) )*(DEXP(X(2)+X(4) )+ 

+  ( 1/DEXP(X(2)+X(4) ) )+(DEXP(X(2) )/ 

+  DEXP(X(4) ) )+(DEXP(X(4) )/DEXP(X(2) ) ) )-LOG(23 .00) 

*******************. '.■********■/.**  ************************ 


****** 


4  ST  EQUATION 

F(4)=DEXP(X(1)+X(4))*(DEXP(X(2)+X(3))+ 

+  ( 1/DEXP(X(2)+X(3) ) )+(DEXP(X(2) )/ 

+  DEXP(X(3) ) )+(DEXP(X( 3) )/DEXP(X(2) ) ) )-LOG( 14 .00) 

RETURN 
END 
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APPENDIX  F.  PROGRAM  LLMDEP 


PROGRAM  LLMDEP 


THIS  IS  THE  FORTRAN  PROGRAM  TO  CALCULATE  THE  RELIABILITY 
OF  PYROTECHNIC  DEVICE.  IT  ASSUMES  THAT  THERE  IS  A  DEPENDENCE 
BETWEEN  MANUFACTURER  AND  ENVIRONMENT  TESTS.  EXPECTATION-MAXIMIZA. 
TION  ALGORITHM  IS  USED  IN  THIS  MODULE.  ALGORITHM  STARTS  WITH 
INITIAL  GUESSES  FOR  PARAMETERS  AND  ESTIMATES  EXPECTATIONS.  IT 
RECALCULATES  CELL  PROBABILITIES  AND  UPDATES  EXPECTATIONS  UNTIL 
IT  CONVERGES.  AN  ITERATIVE  NEWTON  AND  RAPHSON  PROCEDURE  IS  USED 
DURING  UPDATATION  OF  CELL  PROBABILITIES.  THIS  PROCEDURE  IS  DONE 
BY  A  SUBROUTINE  NAMED  UCPROB . 


VARIABLES 


FOM  :  NUMBER  OF  FAILURES  IN  MANUFACTURER  TEST. 

FOTH  :  NUMBER  OF  FAILURES  IN  TEMP.  AND  HUMIDITY  TEST. 

FOV  :  NUMBER  OF  FAILURES  IN  VIBRATION  TEST. 

FOA  :  NUMBER  OF  FAILURES  IN  ALTITUDE  TEST. 

MU  :  OVERALL  MEAN. 

LP1  :  MEAN  EFFECT  OF  MANUFACTURER  TEST. 

LP2  :  MEAN  EFFECT  OF  MANUFACTURER  TEST. 

LP3  :  MEAN  EFFECT  OF  MANUFACTURER  TEST. 

LP4  :  MEAN  EFFECT  OF  MANUFACTURER  TEST. 

TETHA  :  TWO  WAY  INTERACTION  TERMS 
RHMLE  :  RELIABILITY  OF  DEVICE 

MPPPP  :  CELL  FREQUENCY  WITH  RESPECT  TO  TESTS  RESULTS. 

SUM  OF  CELL  FREQUENCIES  WHICH  HAVE  PASSED  DEVICES 
MANUFACTURER  AND  VIBRATION  TEST. 


MPOPO 


MPOOO  :  SUM  OF  CELL  FREQUENCIES  WHICH  HAVE  PASSED  DEVICES 


*  MANUFACTURER  TEST. 

*  IP  :  INITIAL  PROBABILITY  VECTOR 

*  FP  :  UPDATED  (FINAL)  PROBABILITY  VECTOR 

*  Y  :  CELL  EXPECTATION  VECTOR 

*  A,B,C  :  SOME  TERMS  TO  MAKE  THE  CALCULATIONS  EASY. 

*  EXPPPP  :  EXPECTED  NUMBER  OF  DEVICES  IN  CELL  WHICH 

*  HAS  A  RESULTANT  VECTOR  (  P  P  P  P  )  IN  MANUFACTURER, 

*  TEMPERATURE,  VIBRATION  AND  ALTITUDE  TEST  RESPECTIVELY 

*  FLAG  :  INDICATOR  VARIABLE  OF  CONVERGENCE  FOR  PARAMETERS. 

*  RHMLE  :  ESTIMATED  RELIABILITY  OF  PYROTECHNIC  DEVICE. 

*  TYPE  DECLERATION 
PARAMETER  (K=10000) 

LOGICAL  FLAG (16) 

INTEGER  I 

REAL  FOM,  FOTH ,  FOV,  FOA,  LP1,  LP2 ,  LP3,  LP4,  TETHA 
REAL  MPPPP , MPPPF , MPPFP , MPPFF , 

+  MPFPP , MPFPF , MPFFP , MPFFF , 

+  MFPPP , MFPPF , MFPFP , MFPFF , 

+  MFFPP , MFFPF , MFFFP , MFFFF 
REAL  EXPPPP(K) ,EXPPPF(K) .EXPPFP(K) ,EXPPFF(K) , 

+  EXPFPP(K) .EXPFPF(K) .EXPFFP(K) ,EXPFFF(K) , 

+  EXFPPP(K) .EXFPPF(K) ,EXFPFP(K) ,EXFPFF(K) , 

+  EXFFPP ( K ) , EXFFPF ( K ) , EXFFFP (K) , EXFFFF ( K ) 

REAL  MPOOO,  MFOOO ,  MPPOO,  MPOPO,  MPOOP,  A,  B,  C 
REAL  IP( 16 , 1) ,FP( 16,1),X(16) 

REAL  Y( 16) 

COMMON  /  PROB  /  Y 

CALL  EXCMS  (’FILEDEF  IS  DISK  END  DATA  A1  '  ) 

*  INITIALIZATION 
RHMLE  =0.0 
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10 


EPS  =  0.001 
DO  10  L=1 , 16 

FLAG(L)=. FALSE. 
CONTINUE 


READING  THE  NUMBER  OF  FAILURES  IN  EACH  TEST  INTERACTIVELY 
WRITE(*,*) 'PLEASE  ENTER  THE  #  OF  FAILURES  IN  MANUFACTURER  TEST’ 
READ(* ,*)  FOM 

WRITE (*,*)’ PLEASE  ENTER  THE  it  OF  FAILURES  IN  TEMP.  AND  HUM.  TEST' 
READC*,*)  FOTH 

WRITE(*,*) 'PLEASE  ENTER  THE  it  OF  FAILURES  IN  VIBRATION  TEST' 
READ(* ,*)  FOV 

WRITE(*,*) 'PLEASE  ENTER  THE  it  OF  FAILURES  IN  ALTITUDE  TEST' 
READ(*,*)  FOA 


READING  THE  INITIAL  GUESS  FOR  EACH  CELL  IN  HYPOTHETICAL 
CONTINGENCY  TABLE 

WRITE  (* ,*) ' PLEASE  ENTER  INITIAL  GUESS  FOR  MPPPP’ 
READ(*,*)  MPPPP 

WRITE  (*,*)' PLEASE  ENTER  INITIAL  GUESS  FOR  MPPPF’ 

READ(* ,*)  MPPPF 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MPPFP' 

READ(* ,*)  MPPFP 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MPPFF' 

READ(* ,*)  MPPFF 

WRITE  (*,*)' PLEASE  ENTER  INITIAL  GUESS  FOR  MPFPP' 

READ(* ,*)  MPFPP 

WRITE  (*,*)' PLEASE  ENTER  INITIAL  GUESS  FOR  MPFPF' 

READ(* ,*)  MPF^F 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MPFFP ' 
READ(*,*)  MPFFP 

WRITE  (*,*)' PLEASE  ENTER  INITIAL  GUESS  FOR  MPFFF ' 
READ(*,*)  MPFFF 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFPPP' 
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READ(*,*)  MFPPP 

WRITE  (*,*)' PLEASE  ENTER  INITIAL  GUESS  FOR  MFPPF* 

READ(* ,*)  MFPPF 

WRITE  (*,*)' PLEASE  ENTER  INITIAL  GUESS  FOR  MFPFP * 

READ(* ,*)  MFPFP 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFPFF ' 

READ(* ,*)  MFPFF 

WRITE  (*,*) ' PLEASE  ENTER  INITIAL  GUESS  FOR  MFFPP' 

READ(* ,*)  MFFPP 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFFPF ' 

READ(* ,*)  MFFPF 

WRITE  (*,*)' PLEASE  ENTER  INITIAL  GUESS  FOR  MFFFP1 
READ(* ,*)  MFFFP 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFFFF ’ 

READ(*,*)  MFFFF 

**********  A-frb*irk**1rM*tt**±**it**1rt.*+.***********i:**±****-ic-k-i:** 

CALCULATION  OF  INITIAL  PROBABALITIES  USING  CELL  FREQUENCIES 

IP(1,1)  =  MPPPP/92 . 00 

IP(2,1)  =  MPPPF/92 . 00 

IP( 3,1)  =  MPPFP/92 . 00 

IP(4, 1)  =  MPPFF/92 . 00 

IP(5 , 1)  =  MPFPP/92 . 00 

IP(6 , 1)  =  MPFPF/92 . 00 

IP( 7,1)  =  MPFFP/92 . 00 

I P ( 8,1)  =  MPFFF/92 . 00 

IP(9 , 1)  =  MFPPP/92 . 00 

IP(10,1)=  MFPPF/92 . 00 

IP( 11,1)=  MFPFP/92 . 00 

I P ( 12,1)=  MFPFF/92 . 00 

I P ( 13,1)=  MFFPP/92 . 00 

IP( 14 , 1 )=  MFFPF/92 . 00 

IP( 15,1)=  MFFFP/92 . 00 
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IP( 16,1)=  MFFFF/92.00 


*  UPDATATION  OF  CELL  FREQUENCIES 
MPPPP=FP( 1 , 1)*92 . 00 
MPPPF=FP(2, 1)*92.00 
MPPFP=FP( 3 , 1 )*92 . 00 
MPPFF=FP(4, 1)*92 . 00 
MPFPP=FP(5,1)*92.00 
M?FFF=FP(6 , 1)*92 . 00 
MPFFP=FP( 7 , 1 )*92 . 00 
MPFFF=FP(8 , 1)*92 . 00 
MFPPP=FP(9 , 1)*92 .00 
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MFPPF=FP( 10 , I)*92 . 00 
MFPFP=FP( 11 , 1)*92 . 00 
MFPFF=FP( 12 , 1)*92 . 00 
MFFPP=FP( 13 , 1)*92 . 00 
MFFPF=FP( 14,1)*92.00 
MFFFP=FP( 15 , 1)*92 . 00 
MFFFF=FP( 16 , 1)*92 . 00 


r*** 


INITIAL  EXPECTATIONS 
EXPPPP( 1 )=MPPPP 
EXPPPF( 1)=MPPPF 
EXPPFP( 1 )=MPPFP 
EXPPFF( 1)=MPPFF 
EXPFPP( 1)=MPFPP 
EXPFPF(1)=MPFPF 
EXPFFP( 1)=MPFFP 
EXPFFF( 1)=MPFFF 
EXFPPP(1)=MFPPP 
EXFPPF ( 1 ) =MFPPF 
EXFPFP ( 1 ) =MFPFP 
EXFPFF ( 1 ) =MFPFF 
EXFFPP( 1)=MFFPP 
EXFFPF ( 1 ) =MFFPF 
EXFFFP ( 1 ) =MFFFP 
EXFFFF ( 1 ) =MFFFF 

******  A  A'  k  A  A  *****  ft**  *  ft  *  *  ft  ft  A  ft  A  A  A  A  *********  A  A  ft-A  Aftft  A-**A  A  Aft  A  Aft-ft  ****** 

MPOOO=MPPPP+MPPPF+MPPFP+MPPFF+MPFPP+MPFPF+MPFFP+MPFFF 

MFOOO=MFPPP+MFPPF+MFPFP+MFPFF+MFFPP+MFFPF+MFFFP+MFFFF 

MPPOO=MPPPP+MPPPF+MPPFP+MPPFF 

MPOPO=MPPPP+MPPPF+MPFPP+MPFPF 

MPOOP=MPPPP+MPPFP+MPFPP+MPFFP 

A=92 . 00-MPP00 

B=92 . OO-MPOPO 
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C=92 . OO-MPOOP 


NEXT  EXPECTATIONS 
DO  20  1=2, K 

EXPPPP( I )=( 20 . 00-FOM)*(MPPPP/MP000 )+( 20 . 00-FOTH)*(MPPPP/MPP00 )+ 
+  ( 32 . 00-FOV)*(MPPPP/MP0P0)+( 20 . 00-FOA)*(MPPPP/MPOOP) 

IF(ABS(EXPPPP(I)-EXPPPP(I-1)) .LE.EPS)  THEN 
FLAG( 1)=.TRUE. 

END  IF 

EXPPPF(I)=(20.00-FOM)*(MPPPF/HP000)+(20.00-FOTH)*(MPPPF/MPP00)+ 
+  ( 32 . 00-FOV)*(MPPPF /MPOPO )+FOA*(MPPPF/C) 

IF( ABS(EXPPPF( I ) -EXPPPF( I- 1 ) ) . LE.EPS)  THEN 
FLAG(2)=.TRUE. 

END  IF 

EXPPFP(I)=(20.00-FOM)*(MPPFP/MP000)+C20.00-FOTH)*(MPPFP/MPP00)+ 
+  FOV*(MPPFP/B)+( 20 . 00-FOA)*(MPPFP/MP00P) 
IF(ABS(EXPPFP(I)-EXPPFP(I-1)). LE.EPS)  THEN 
FLAG(3)=.TRUE. 

END  IF 

EXPPFF(I)=(20.00-FOM)*(MPPFF/MP000)+(20.00-FOTH)*(MPPFF/MPP00)+ 
+  FOV* ( MPPFF / B ) +FOA* ( MPPFF /C ) 

I F ( ABS ( EXPPFF ( I ) - EXPPFF (1-1)). LE.EPS)  THEN 
FLAG(4)= . TRUE . 

END  IF 

EXPFPP( I )=( 20 . 00-FOM)*(MPFPP/MP000)+FOTH*(MPFPP/A)+( 32 . OO-FOV)* 
+  (MPFPP/MPOP0)+(2O.00-FOA)*(MPFPP/MP00P) 

IF( ABS(EXPFPP( I )-EXPFPP( I- 1) ) .LE.EPS)  THEN 
FLAG(5)=.TRUE. 

END  IF 

EXPFPF( I )=( 20 . 00-FOM)*(MPFPF/MP000 )+FOTH*(MPFPF/A)+( 32 . OO-FOV)* 
+  ( MPFPF/MPOPO ) +FOA*( MPFPF/C ) 

IF(ABS(EXPFPF(I)-EXPFPF(I-1)). LE.EPS)  THEN 
FLAG(6)=.TRUE. 

END  IF 
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EXPFFP( I )=(20 . 00-FOM)*(MPFFP/MP000 )+FOTH*(MPFFP/A)+FOV*(MPFFP/B) 
( 20 . 00-FOA)*(MPFFP/MP00P) 

IF( ABS(EXPFFP( I ) -EXPFFP( I -  1 ) ) . LE . EPS)  THEN 
FLAG( 7 )= . TRUE . 

END  IF 

EXPFFF(I)=(20 ,00-FOM)*(MPFFF/MPOOO)+FOTH*(MPFFF/A)+FOV*(MPFFF/B) 
FOA*(MPFFF/C) 

IF(ABS(EXPFFF( I ) -EXPFFF( I - 1 ) ) . LE . EPS )  THEN 
FLAG( 8)= . TRUE . 

END  IF 

EXFPPP( I ) =FOM*(MFPPP/MFOOO ) +FOTH*( MFPPP/A ) +FOV*( MFPPP/B )+FOA* 
(MFPPP/C) 

IF(ABS(EXFPPP( I ) -EXFPPP( I- 1) ) . LE .EPS)  THEN 
FLAG(9)=.TRUE. 

END  IF 

EXFPPF ( I ) =FOM* ( MFP  PF/ MFO  0  0 ) +FOTH* (MFPPF/A) +FOV* (MFPPF/B) +FOA* 
(MFPPF/C) 

IF(ABS(EXFPPF( I ) -EXFPPF( I- 1) ) . LE .EPS)  THEN 
FLAG( 10)=. TRUE. 

END  IF 

EXFPFP ( I ) =FOM* (MFPFP/MF000) +FOTH* ( MFPFP / A ) +FOV* (MFPFP/B) +FOA* 
(MFPFP/C) 

IF(ABS(EXFPFP( I ) -EXFPFP(I-l) ) . LE . EPS )  THEN 
FLAG( 11)=. TRUE. 

END  IF 

EXFPFF ( I ) =FOM* ( MFPFF/MFO 0  0 ) +FOTH* ( MFPFF/A ) +FOV* (MFPFF/B) +FOA* 
(MFPFF/C) 

IF(ABS(EXFPFF( I ) - EXFPFF ( I - 1 ) ) . LE . EPS )  THEN 
FLAG( 12)=. TRUE. 

END  IF 

EXFFPP( I ) =FOM*( MFFPP/MFOO  0 ) +FOTH* ( MFFPP/A ) +FOV*( MFFPP/ B ) +FOA* 
(MFFPP/C) 

IF(ABS(EXFFPP(I)-EXFFPP(I-1)) .LE.EPS)  THEN 


FLAG( 13)=. TRUE. 


END  IF 

EXFFPF ( I ) =FOM* ( MFFPF/ MFO  0  0 ) +FOTH* ( MFFPF/A ) +FO V - ( MFFPF / B ) +FOA* 
+  (MFFPF /C) 

IF(ABS(EXFFPF( I ) -EXFFPFC I- 1 ) ) . LE . EPS)  THEN 
FLAG( 14)= . TRUE . 

END  IF 

EXFFFP(I)=FOM*(MFFFP/MFOOO)+FOTH*(MFFFP/A)+FOV*(MFFFP/B)+FOA* 
+  (MFFFP/C) 

I F ( ABS ( EXFFFP ( I ) - EXFFFP ( I-1)).LE.EPS)  THEN 
FLAG( 15 )= . TRUE . 

END  IF 

EXFFFF ( I ) =FOM*( MFFFF/MFO  0 0 ) +FOTH* ( MFFFF/A) +FOV* ( MFFFF/B ) +FOA* 
+  (MFFFF/C) 

I F ( ABS ( EXFFFF ( I ) - EXFFFF (I-1)).LE.EPS)  THEN 
FLAG ( 1 6 ) = . TRUE . 

END  IF 

MPPPP  =  EXPPPP(I) 

MPPPF  =  EXPPPF(I) 

MPPFP  =  EXPPFP(I) 

MPPFF  =  EXPPFF(I) 

MPFPP  =  EXPFPP(I) 

MPFPF  =  EXPFPF(I) 

MPFFP  =  EXPFFP(I) 

MPFFF  =  EXPFFF ( I ) 

MFPPP  =  EXFPPP(I) 

MFPPF  =  EXFPPF(I) 

MFPFP  =  EXFPFP(I) 

MFPFF  =  EXFPFF(I) 

MFFPP  =  EXFFPP(I) 

MFFPF  =  EXFFPFC I) 

MFFFP  =  EXFFFP(I) 
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MFFFF  =  EXFFFF(I) 


CHECK  FOR  THE  STOPING  CONDITION 

IF(  FLAG(l) .AND.FLAC(2) . AND.FLAG(3) . AND.FLAG(4) . AND.FLAG(5) 

+  . AND . FLAG ( 6 ) . AND . FLAG ( 7 ) . AND . FLAG ( 8 ) . AND .FLAG( 9 ) . AND . FLAG( 10) 
+  . AND . FLAG ( 1 1 ) . AND . FLAG ( 1 2 ) . AND . FLAG ( 1 3 ) . AND . FLAG ( 14 ) . AND . 

+  FLAG (15). AND. FLAG (16))  THEN 


CALCULATION  OF  FINAL  EXPPPP  (STOPING  CONDITION  IS  MET) 
RHMLE=EXPPPP( I ) /92 . 00 
GO  TO  30 
END  IF 


STOPING  CONDITION  IS  NOT  MET.  PROBABILITIES  FOR  THE  NEXT 
NEWTON  AND  RAPHSON  PROCEDURE 
IP(1,1)  =  MPPPP/92 . 00 
IP(2 , 1)  =  MPPPF/92 . 00 
IP( 3,1)  =  MPPFP/92 . 00 
IP(4 , 1)  =  MPPFF/92.00 
IP(5 , 1)  =  MPFPP/92 . 00 
IP(6 , 1)  =  MPFPF/92 . 00 
IP( 7,1)  =  MPFFP/92 . 00 
IP(8, 1)  =  MPFFF/92 . 00 
IP(9 , 1)  =  MFPPP/92.00 
IP(10,1)=  MFPPF/92  00 
IP(11,1)=  MFPFP/92 . 00 
IP( 12,1)=  MFPFF/92 . 00 
I P ( 13,1)=  MFFPP/92 . 00 
I P  ( 14,1)=  MFFPF/92 . 00 
IP( 15,1)=  MFFFP/92 . 00 
IP( 16,1)=  MFFFF/92 . 00 


CALL  UCPROB(IP.FP) 
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MPPPP=FP( 1 , 1 )*92 . 00 
MPPPF=FP( 2 , 1 )*92 . 00 
MPPFP=FP( 3 , 1)*92 . 00 
MPPFF=FP(4 , 1 )*92 . 00 
MPFPP=FP(5 ,  1)*92 . 00 
MPFPF=FP(6, 1)*92.00 
MPFFP=FP( 7 , 1 )*92 . 00 
MPFFF=FP( 8,1)*92.00 
MFPPP=FP( 9 , 1)*92 . 00 
MFPPF=FP( 10 , 1)*92 . 00 
MFPFP-F P( 11 , 1)*92 . 00 
MFPFF=FP( 12,1)*92.00 
MFFPP=FP( 13, 1)*92 .00 
MFFPF=FP(14,1)*32.00 
MFFFP=FP( 15 , 1)*92 . 00 
MFFFF=FP( 16 , 1 )*92 . 00 

MPOOO=MPPPP+MPPPF+MPPFP+MPPFF+MPFPP+MPFPF+MPFFP+MPFFF 

MFOOO=MFPPP+MFPPF+MFPFP+MFPFF+MFFPP+MFFPF+MFFFP+MFFFF 

MPPOO=MPPPP+MPPPF+MPPFP+MPPFF 

MPOPO=MPPPP+MPPPF+MPFPP+MPFPF 

MPOOP=MPPPP+MPPFP+MPFPP+MPFFP 

A=92 . 00-MPP00 

B=92 . OO-MPOPO 

C=92 . OO-MPOOP 

20  CONTINUE 

30  WR ITE ( 15 , 40 )FOM , FOTH , FOV , FOA , RHMLE 

40  FORMAT ( / , 5X , ' CASE ' ,4X ,F5 . 2 , 2X ,F5 . 2 , 2X ,F5 . 2 , 2X ,F5 . 2 , / , 15X , 
+'MLE  =  ' ,F12 . 7 ) 


STOP 
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END 


SUBROUTINE  UCPROB  (IP,FP) 


THIS  SUBROUTINE  UPDATES  CELL  PROBABILITIES  USING  NEWTON  AND 


RAPHSON  PROCEDURE  WHICH  IS  DESCRIBED  IN  SAS  . 


*  VARIABLES 

*  FOM  :  NUMBER  OF  FAILURES  IN  MANUFACTURER  TEST. 

*  FOTH  :  NUMBER  OF  FAILURES  IN  TEMP.  AND  HUMIDITY  TEST. 

*  FOV  :  NUMBER  OF  FAILURES  IN  VIBRATION  TEST. 

*  FOA  :  NUMBER  OF  FAILURES  IN  ALTITUDE  TEST. 

*  MU  :  OVERALL  MEAN. 

*  LP1  :  MEAN  EFFECT  OF  MANUFACTURER  TEST. 

*  LP2  :  MEAN  EFFECT  OF  MANUFACTURER  TEST. 

*  LP3  :  MEAN  EFFECT  OF  MANUFACTURER  TEST. 

*  LP4  :  MEAN  EFFECT  OF  MANUFACTURER  TEST. 

*  TETHA  :  TWO  WAY  INTERACTION  TERMS 

*  RHMLE  :  RELAIBILITY  OF  DEVICE 

*  MPPPP  :  CELL  FREQUENCY  WITH  RESPECT  TO  TESTS  RESULTS. 

*  MPOPO  :  SUM  OF  CELL  FREQUENCIES  WHICH  HAVE  PASSED  DEVICES 

*  MANUFACTURER  AND  VIBRATION  TEST. 

*  MPOOO  :  SUM  OF  CELL  FREQUENCIES  WHICH  HAVE  PASSED  DEVICES 

*  MANUFACTURER  TEST. 

*  IP  :  INITIAL  PROBABILITY  VECTOR 

*  FP  :  UPDATED  (FINAL)  PROBABILITY  VECTOR 

*  Y  :  CELL  EXPECTATION  VECTOR 

*  A  ,B  ,C  :  SOME  TERMS  TO  MAKE  THE  CALCULATIONS  EASY. 

*  EXPPPP  :  EXPECTED  NUMBER  OF  DEVICES  IN  CELL  WHICH 

*  HAS  A  RESULTANT  VECTOR  (  P  P  P  P  )  IN  MANUFACTURER, 

*  TEMPERATURE,  VIBRATION  AND  ALTITUDE  TEST  RESPECTIVELY 

*  SIGN  :  INDICATOR  VARIABLE  OF  CONVERGENCE  FOP  PARAMETERS. 

*  RHMLE  :  ESTIMATED  RELIABILITY  OF  PYROTECHNIC  DEVICE. 
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•***: 


TYPE  DECLERATION 
LOGICAL  SIGN(5) 

REAL  IP( 16 , 1) ,FP( 16,1) 

REAL  F0(15,1),F1(15,1),S(15,15),X(15,5),B0(5,1), 
+B1(5,1),PI0(15,1),PI1(15,1) ,C(5,5) ,G(5,1) ,SINV (15,15) ,F( 15 , 1) , 
+PR|(15,1),PR2(5,1),PR3(15,5),PR4(5,5), SUM .LAST, PIO 1( 16,1) ,LHE , 
+LHEMAX ,RIP( 15,1), PR 5(15, 5) ,CINV(5,5) ,DIF( 15 , 1 ), PR6(5 , 1) .LAMBDA , 
+PI 11(16,1) , EPS, BLAST(5 , 1) ,PLAST( 15,1) ,BNEW(5, 1) ,FLAST( 15, 1) 

+ ,U ,XT(5 , 15) ,PR4INV(5 ,5) 

REAL  YC 16) 

COMMON  /  PROB  /  Y 
INTEGER  I , J ,K 


INITIALIZATION 


DO  50  1=1,5 

SIGN( I )=  .FALSE. 

50  CONTINUE 

DO  70  1=1,5 
DO  60  J=1 , 5 
S(I,J)=0.0 
60  CONTINUE 
70  CONTINUE 

*  *************  A  A  A  *  /;  *  A-  A'  *  *  *  *  'A  ■/.  *  'A  A  Mck*  *  k  **************** 

*  READING  THE  DESIGN  MATRIX 

CALL  EXCMS  ( 'FILEDEF  9  DISK  DESIGN  INPUT  Al') 

DO  80  1=1,15 

READ(9 ,*)  X(I,1),X(I,2),X(I,3),X(I,4),X(I,5) 

80  CONTINUE 
REWIND  9 

*  i?'kmk-kmk+;+:-k  "fir  J? -kit  -A-  A'  A'  f/A-A-  A-  A  A  A*  A*  A  A  A~A  A  /.rAA'A-AA  A  A' A  A  A  A  A  A  AAA A A AAA AA  A  AA’ 

*  INVERSE  OF  VARIANCE  AND  COVERIANCE  MATRIX  FOR  INITIAL  BO 
DO  100  1=1,15 
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DO  90  J=1 , 15 

IF(I.EQ.J)  THEN 

SINV( I ,I)=(IP(I,1)-(IP(I, 1)**2 . 0) )*92 . 00 
END  IF 

IF(I.NE.J)  THEN 

SINV(I,J)=(-IP(J,1))*IP(I,1)*92.00 
END  IF 

90  CONTINUE 
100  CONTINUE 


LOGIT  RESPONSE  FUNCTIONS 
DO  110  1=1,15 

F(I , 1)  =AL0G(IP(I,1)/IP(16,1)) 
110  CONTINUE 


*  THE  TRANSPOZE  OF  THE  DESIGN  MATRIX 
CALL  TRNRR  ( 15 ,5 ,X, 15 ,5 , 15 ,XT,5) 

*  ftftft-AA  A  ft  Aft  A  ft  Aftftftftftftft^ftftft  ft  ftftftft  ft  A  ftftftft ftftftftftft  ft  ftftftftftftftftfcftftft-ftftftftftftftft-ft-ftft 

*  MATRIX  MULTIPLICATION  PR1=(SINV*F) 

CALL  MRRRR  ( 15 , 15 ,SINV , 15 , 15 , 1 ,F, 15 , 15 , 1 ,PR1 , 15) 

ft  ft  ft*  ft  A  A  ft  ft  ft  A  A  ftftftftftftftftftftftft  ft  Aft  A  A  A  A  ftftftftftftftftftftftftftft-ftftft  ft  ftftftftftft  ftftftftft  ftftftft  A 

*  MATRIX  MULTIPLICATION  PR2=(XT*PR1) 

CALL  MRRRR  (5 , 15 ,XT,5 , 15 , 1,PR1 , 15 ,5 , 1 ,PR2 ,5) 

ft  ftftftftftft*  ft  ft  A  A  ft ftftftftft-ftftftftftftft  ftftft  ftft  ft  ft-ftftftftft*ftftftftftftftftftftftftftft*ftftftftft*ftftftftftft 

*  MATRIX  MULTIPLICATION  PR3=(SINV*X) 

CALL  MRRRR  ( 15 , 15 ,SINV, 15 , 15 ,5 ,X , 15 , 15 ,5 ,PR3, 15) 

ft  ft  ft-ftftftftft  ftftftft  ft  A  ftftft  ftftftftftft  A  A  ftftftftft1  *■  ftftftft*  Aft1  ftft* *******  ft  ftftftftft  ftftftftft  ft  ft***-* 

*  MATRIX  MULTIPLICATION  PR4=(XT*PR3) 

CALL  MRRRR  (5 , 15 ,XT,5 , 15 ,5 ,PR3 , 15 ,5 ,5 ,PR4 ,5) 

ft  ft  ftftftftft  ft  A  ft  ft  ft  ft  ft  ft  Aft  ftftft  ft  ftftftft'ftftft  ft  ft  A  ftftftft  ft  ft  ftft  ft  ft  ftftftftftft  ft  ftftft-ftftftft  ftftftft  A  ftftftft  ft 

*  INVERSE  OF  THE  MATRIX  MULTIPLICATION  PR4=PR4INV 

CALL  LINRG  (5,PR4,5,PR4INV,5) 

ft  ftftftftftft*  ftft  ft  ft  ftftftftftftftft  ft-ft-ft-ftft  ftft-ftft  ft  ft  **********  ftftftftftft*  ft  ft  A-A-ftft'ftftftftftft'A  ftftft* 
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INITIAL  ESTIMATION  OF  PARAMETERS  BO 
MATRIX  MULTIPLICATION  B0=(PR4INV*PR2) 

CALL  MRRRR  (5,5, PR4INV ,5,5,1, PR2 , 5 , 5 , 1 , BO , 5 ) 


FO=  X  *  BO 

CALL  MRRRR  ( 15 , 5 ,X , 15 , 5 , 1 ,BO , 5 , 15 , 1 ,FO , 15) 


INITIAL  PROBABILITIES  PI0=(EXP(F0 ) ) 
DO  120  1=1,15 

PIO( I , l)=EXP(FO(I , 1) ) 

120  CONTINUE 
SUM=0 . 0 
DO  130  1=1,15 

SUM=SUM+PIO( I , 1) 

130  CONTINUE 

LAST=1 .0/(1. O+SUM) 


PROBABILITY  MATRIX  WHICH  INCLUDES  16  VALUES  PIOl 
FOR  THE  INITIAL  ESTIMATE  OF  LIKELIHOOD  ESTIMATION 
DO  140  1=1,15 

PI01(I,1)=PI0(I, 1)*LAST 
140  CONTINUE 

PI01( 16, 1)=LAST 


INITIAL  LIKELIHOOD  FOR  NEXT  ITERATION  AT  STEP  BO 
LHE=0 . 0 
DO  150  1=1,16 

LHE=LHE+Y(I )*ALOG(PI01( I , 1) ) 

150  CONTINUE 
LHEMAX=LHE 


REORGANIZED  INITIAL  PROBABILITIES  FOR  UPDATATION  RIP 
DO  160  1=1,15 

RIP( I , 1 )=IP( 1,1) 
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160  CONTINUE 


fcL  I kJUJL  JLJtwJLJLp  JUJLJU > I ^  «JU  ijL  ■!■  < U «^U X ai < J*  J** I*  J * «J »  Jj  #V"<VvV  '/»  ^w^h|U ^  X  JLJ^k^i  jU X 

*  FIRST  ITERATION  IN  NEWTON  AND  RAPHSON  METHOD 
DO  170  1=1,15 

DO  180  J=1 , 15 

IF(I.EQ.J)  THEN 

SINV(  I  ,I)=(PI0(I,1)-(PI0(1, 1)**2 . 0) )*92 . 00 
END  IF 

IF(I.NE.J)  THEN 

SINV(I,J)=(-PI0(J,1))*PI0(I,1)*92.00 
END  IF 

170  CONTINUE 
180  CONTINUE 

*  tt-**±-irirtt*f.*iri;**rfrfi!rK*ie*1-}r******tt***tt*;rJ.-/rt;**->rk*****±i:*±*-^±****: 

MATRIX  MULTIPLICATION  PR5=(SINV*X) 

CALL  MRRRR  ( 15 , 15 ,SINV, 15 , 15 ,5 ,X , 15 , 15 ,5 ,PR5 , 15) 

MATRIX  MULTIPLICATION  C=(XT*PR5) 

CALL  MRRRR  (5, 15, XT, 5,15, 5, PR5, 15, 5, 5, C, 5) 

**  /.  ft-A  ******  A  **** A  A  *****  A  *****  A  *****  A-  *****  A  *  A  A  A  A  A  A  A  A  A  A  A  A  A  AAA  A*A* 

INVERSE  OF  THE  MATRIX  C=CINV 
CALL  LI NRG  (5 ,C ,5 ,CINV,5) 

*  ******  *****  A  *★**  A  A  **********  A  A  ****  A  A  A'A*  AAA  *****  A******  AAA  AAA  **** 

DO  190  1=1,15 

DIF( I , 1)=  92  00*(RIP(I,1)-PI0(I,1)) 

190  CONTINUE 

*  **  ft  ft**  ***********  A  ******  *****  A  *****  A  *■*■*  *  A r***********r************ 

*  MATRIX  MULTIPLICATION  FOR  G=(XT*DIF) 

CALL  MRRRR  (5 , 15 ,XT,5 , 15 , 1 ,DIF, 15 ,5 , 1 ,G,5) 

*  ******  A*  AAA  *  **A*i  A  ******* '.'."A  A  A  A  A  A  "A1  A  AAA  **■*"*  *****-**  *  A  A A  A  A  A  A  A  A'  A  A  *■*■*■* 

*  MATRIX  MULTIPLICATION  (PR6=  CINV*G) 

CALL  MRRRR  (5 ,5 ,CINV ,5 ,5 , 1 ,G ,5 ,5 , 1 ,PR6 ,5  ) 

LAMBDA=1 .0 


* 


* 

* 


* 
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200  DO  210  1=1,5 

PR6 ( I , 1 ) =PR6 ( I , 1)* *LAMBDA 
210  CONTINUE 


INITIAL  VALUE  FOR  B1 
DO  220  1=1,5 

!( I , 1 )=B0( I , 1) -PR6( 1,1) 
220  CONTINUE 


MATRIX  MULTIPLICATION  F1=(X*B1) 

CALL  MRRRR  ( 15 ,5 ,X, 15 ,5 , 1 ,B1 ,5 , 15 , 1 ,F1 , 15) 


CALCULATION  OF  PROBABILITIES  FOR  B1  PI1=EFP(F1) 
DO  230  1=1,15 

PI 1 ( 1 , 1 )=EXP(F1( 1 , 1) ) 

230  CONTINUE 


*  CALCULATION  OF  THE  16  TH  PROBABILITY  VALUE 
SUM=0 . 0 

DO  240  1=1,15 

SUM=SUM+PI 1(1,1) 

240  CONTINUE 

LAST=1 .0/(1. 0+SUM) 

*  ************************** 


DO  250  1=1,15 

PI11(I,1)=PI1(I,1)*LAST 
250  CONTINUE 

PI  11 ( 16, 1)=LAST 


*  INITIAL  LIKELIHOOD  ESTIMATION 
LHE=0 . 0 
DO  260  1=1,16 

LHE=LHE+Y ( I )*ALOG(PI 11(1,1)) 
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260  CONTINUE 


EPS=0 . 00 1 

I F ( LHE . LT . LHEMAX )  THEN 
K=K+1 

IF(K.GT.IO)  THEN 
DO  270  1=1,5 

BLAST ( 1 , 1 ) =B0 ( 1 , 1 ) 
CONTINUE 
DO  280  1=1,15 

PLAST( I, 1)=P 101(1,1) 
CONTINUE 
GO  TO  490 


END  IF 

LAMBDA=LAMBDA/2.0 
GO  TO  200 


END  IF 
LHEMAX=LHE 


DO  290  1=1,5 


IF(ABS(B0( I ,1)-B1(I,1)) .LE .EPS)  THEN 
SIGN(I)=.TRUE. 

END  IF 


290  CONTINUE 


CHECKING  CRITERIAS 

IF(SIGN( 1) .AND.SIGN(2) .AND .SIGN(3) .AND.SIGN(4) .AND.SIGN(5))  THEN 
DO  300  1=1,5 

BLAST( I , 1)=B1( I , 1) 

CONTINUE 
DO  310  1=1,15 

PLAST( I , 1)=PI 11(1 , 1) 

CONTINUE 
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GO  TO  490 


END  IF 

CRITERIAS  ARE  NOT  MET  THEN  NEW  ITERATIONS 
DO  320  1=1,5 

BLAST ( I,1)=B1(I,1) 

320  CONTINUE 

DO  330  1=1,15 

PLAST(I,1)=PI1(I,1) 

330  CONTINUE 


340  DO  360  1=1,15 

DO  350  J=1 , 15 

IF(I.EQ.J)  THEN 

SINV(I,I)=(PLAST(I,1)-(PLAST(I,1)**2.0))  *92.00 
END  IF 

IF(I.NE.J)  THEN 

SINV(I,J)=( -PLAST( J , 1) )*PLAST( I , 1 )*92 . 00 
END  IF 

350  CONTINUE 
360  CONTINUE 

*  MATRIX  MULTIPLICATION  PR3=(SINV*X) 

CALL  MRRRR  ( 15 , 15 ,SINV , 15 , 15 ,5 ,X , 15 , 15 ,5 ,PR3 , 15) 


MATRIX  MULTIPLICATION  C=(XT*PR3) 

CALL  MRRRR  (5 , 15 ,XT,5 , 15 ,5 ,PR3 , 15 ,5 ,5 , C,5) 


CALL  LINRG  (  5,C,5,CINV,5  ) 


DO  370  1=1,15 

DIF( I , 1 )=  92 . 00*(RIP( I , 1)-PLAST( 1,1)) 
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370  CONTINUE 


*  MATRIX  MULTIPLICATION  G=(XT*DIF) 

CALL  MRRRR  (5, 15, XT, 5,15, 1, DIF, 15,5,1 ,G, 5) 


MATRIX  MULTIPLICATION  PR6=(CINV*G) 

CALL  MRRRR  (5 , 5  .CINV.5 ,5 , 1 ,G , 5 ,5 , 1 ,PR6 ,5) 


t**; 


1 i l _ i »-  -  » — i .  - 1  -i  -  - 1  — i--i 


LAMBDA=1 . 0 


NEW  PARAMETER  ESTIMATES 
380  DO  390  1=1,5 

BNEW( I , 1 )=BLAST( I , 1 ) - ( LAMBDA*PR6( 1,1)) 
390  CONTINUE 


■r  *********  *********^ 


*  MATRIX  MULTIPLICATION  F1=(X*BNEW) 

CALL  MRRRR  ( 15 ,5 ,X , 15 ,5 , 1 , BNEW ,5 , 15 , 1 ,F1 , 15) 

DO  400  1=1,15 

PI1( I , 1)=EXP(F1( I , 1) ) 

400  CONTINUE 

*  ***  **  **********************  A  K  ******  A  ***  A  A  *********************** 


SUM=0 . 0 
DO  410  1=1,15 

SUM=SUM+PI 1(1,1) 
410  CONTINUE 

LAST=1 -0/(1 .0+SUM) 


DO  420  1=1,15 

PI 11( 1 , 1)=PI1(1 , 1)*LAST 
420  CONTINUE 

PI 1 1( 16 , 1)=LAST 

*  **************************************************************** 


LIKELIHOOD  ESTIMATION 
LHE=0 . 0 
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DO  430  1=1,16 


LHE=LHE +Y( I )*ALOG(PI 11(1,1)) 
430  CONTINUE 


K=0 

EPS=0 .001 

IF( LHE . LT . LHEMAX)  THEN 
K=K+1 

IF(K.GT.IO)  THEN 
GO  TO  490 
END  IF 

LAMBDA=LAMBDA/2 . 0 
GO  TO  380 
END  IF 
LHEMAX  =LHE 


DO  440  1=1,5 

IF(ABS(B1(I , 1 ) -BLAST( 1,1)) . LE . EPS)  THEN 
SIGN( I )= . TRUE . 

END  IF 


440  CONTINUE 


IF(SIGN( 1) .AND . SIGN(2) .AND . SIGN(3) . AND. SIGN(4) .AND . SIGN(5) )  THEN 
DO  450  1=1,5 

BLAST (1,1) =BNEW (1,1) 

CONTINUE 
DO  460  1=1,15 

PLAST( I , 1 ) =P 111(1,1) 

CONTINUE 
GO  TO  490 
END  IF 


CR ITER IAS  ARE  NOT  MET  THEN  NEW  ITERATIONS 
DO  470  1=1,5 
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BLAST( I ,1)=B 1(1,1) 
470  CONTINUE 

DO  480  1=1,15 

PLAST( I, 1)=P 111(1,1) 
480  CONTINUE 
GO  TO  340 


490  CONTINUE 


MATRIX  MULTIPLICATION  FLAST=(X*BLAST) 

CALL  MRRRR  ( 15 , 5 ,X , 15 , 5 , 1 , BLAST ,5 , 15 , 1 ,FLAST, 15) 


CALCULATION  OF  FINAL  PROBABILITIES 
DO  500  1=1,15 

PLAST( I , 1)=EXP(FLAST( 1,1)) 


500  CONTINUE 


SUM=0 . 0 


DO  510  1=1,15 

SUM=SUM+PLAST( 1,1) 


510  CONTINUE 


LAST=1 . 0/( 1 .0+SUM) 

DO  520  1=1,15 

FP( I , 1 ) =PLAST( I , 1 )*LAST 


520  CONTINUE 


FP( 16 , 1)=LAST 
RETURN 


122 


APPENDIX  G.  PROGRAM  MLEB 


PROGRAM  MLEB 


THIS  IS  THE  PROGRAM  TO  CALCULATE  RELIABILITY  OF  THE  DEVICE 
WITH  DEPENDENCE  ASSUMPTION . PROGRAM  ASSUMES  THAT  FAILED  ITEM  FROM 
ANY  OF  ENVIRONMENT  TESTS  FAILS  FROM  MANUFACTURER  TEST  TOO.  THIS 
IS  WORST  CASE  SCENARIO.  IT  READS  NUMBER  OF  SUCCESSFUL  ITEMS  FROM 
AN  INPUT  DATA  CALLED  SUSVECT.  FINALLY  THE  PROGRAM  WRITES  RESULTS 
TO  AN  OUTPUT  FILE  CALLED  RESULT. 


*  TYPE  DECLARATION 

REAL  SOM(5000) ,SOTH(5000) ,SOV(5000) ,SOA(5000) ,X(4) ,R1H ,R2H ,R3H , 
+  R4H ,RHMLE(5000) 
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INTEGER  I 


FILES  FOR  READING  AND  WRITING 

CALL  EXCMS  ( 'FILEDEF  9  DISK  SUSVECT  DATA  A1 ' ) 

CALL  EXCMS  ('FILEDEF  17  DISK  RESULT  DATA  Al') 


READING  NUMBER  OF  SUCCESS  IN  EACH  TEST 
DO  10  1=1,5000 

READ(9 ,*)  SOM( I ) , SOTH( I ) , SOV( I ) , SOA( I ) 


NUMBER  OF  SUCCESSES  IN  EACH  TEST 

X( 1 )=  SOM( I )+SOTH( I )+SOV( I )+SOA( I ) 
X(2)=  SOTH(I) 

X(3)=  SOV(I) 

X(4)  =  SOA(I) 


CALCULATIONS  SUCCESS  PROBABILITIES  IN  EACH  TEST 
R1H=  X(l)/184.0 
R2H=  X(2)/40 .0 
R3H=  X( 3)/64 . 0 
R4H=  X(4)/40 . 0 

RHMLE(I)=  R1H  *  R2H  *  R3H  *  R4H 


WRITING  RESULTS  TO  AN  OUTPUT  FILE  CALLED  RESULT 
WRITE  (17,*)  RHMLE(I) 

CONTINUE 
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APPENDIX  H.  PROGRAM  BONUS 


PROGRAM  BONUS 


THIS  IS  THE  PROGRAM  TO  CALCULATE  BONUS  PERCENTAGE  OF  ANY  FIRM 
WHOSE  LONG  RUN  SUCCESS  PROBABILITIES  ARE  KNOWN.  IN  THIS  PROGRAM 
IT  USES  2000  SUCCESS  VECTORS,  WHICH  ARE  GENERATED  BY  RANVEC  IN 
APPENDIX  B.  THEY  ARE  GENERATED  BY  KNOWN  LONG  RUN  PROBABILITIES 
THE  PROGRAM  USES  TWO  DATA  SETS. THEY  ARE  PRECALCULATED  LCB'S  SETS 
.FIRST  DATA  REPRESENTS  FIRST  INSPECTION,  SECOND  DATA  REPRESENTS 
SECOND  INSPECTION  LOWER  CONFIDENCE  BOUNDS. 

SUCCESS  VECTORS  REPRESENT  OFFERED  LOTS. IT  HAS  A  DETERMINISTIC 
BONUS  LINE. PROGRAM  CALCULATES  LCB  OF  OFFERED  LOT, WITH  FIRST  DATA 
AND  COMPARES  IT  WITH  LCB  OF  BONUS  LINE. IF  FIRM  LCB  IS  GRATER  THAN 
FIRM  GETS  BONUS.  OTHERWISE  FIRM  HAS  A  CHANCE  TO  ONE  MORE  TRY.  IN 
SECOND  TRY,  PROGRAM  CUMULATES  SUCCESS  VECTORS  AND  IT  USES  SECOND 
DATA  TO  FIND  OUT  LCB  OF  CUMULATED  LOT.  AFTER  THIS  CALCULATION  IT 
COMPARES  AGAIN.  FINALLY  IT  COUNTS  NUMBER  OF  TIMES  THAT  THE  FIRM 
GETS  THE  BONUS  IN  1000  REPLICATIONS  AND  ESTIMATES  BONUS  PERCENT. 
IT  WRITES  RESULTS  TO  AN  OUTPUT  FILE  CALLED  BONUS  DATA. 


VARIABLES 


A’A'‘A“A  '/;■  -A-  A'  ’A’  A-  A’A' A-A'A-A-A-ArA  A  A  AAAAAAAAA-AAAAAAAAAAAAAAAAAAAAAAA 


SOM 

SOTH 

SOV 

SOA 

BLINE 

LCB 

FILCB 


NUMBER  OF  SUCCESS  IN  MANUFACTURER  TEST. 

NUMBER  OF  SUCCESS  IN  TEMPERATURE  AND  HUMIDITY  TEST. 
NUMBER  OF  SUCCESS  IN  VIBRATION  TEST. 

NUMBER  OF  SUCCESS  IN  ALTITUDE  TEST. 

LOWER  CONFIDENCE  BOUND  OF  BONUS  LINE 

LOWER  CONFIDENCE  BOUND 

LCB  VALUES  ARRAY  IN  FIRST  INSPECTION 
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* 


•h 


SELCB 

A,B,C,D 

BFI 

BSI 

BTOT 

PRCT 

COUNT 

SIGN 

FLAG 


LCB  VALUES  ARRAY  IN  SECOND  INSPECTION 
DIMENSIONS  FOR  USE  OF  LCB  DATAS 

NUMBER  OF  TIMES  THAT  FIRM  GAT  BONUS  AFTER  1  ST  INSP. 
NUMBER  OF  TIMES  THAT  FIRM  GAT  BONUS  AFTER  2  ST  INSP. 
TOTAL  NUMBER  OF  TIMES  THAT  FIRM  GAT  BONUS. 

BONUS  PERCENT. 

COUNTER  FOR  1000  REPLICATIONS. 

INDICATOR  OF  ACCEPTANCE  FOR  FIRST  INSPECTION 
INDICATOR  OF  ACCEPTANCE  FOR  SECOND  INSPECTION 


TYPE  DECLARATION 
LOGICAL  SIGN(4) , FLAG(4) 

REAL  SOM( 2000 ) , SOTH( 2000 ) ,SOV( 2000) , SOA( 2000) ,BLINE , LCB , BFI , BSI , 
+  BTOT ,FILCB(2 ,2 , 3,2) , SELCB( 3 ,3 ,5 ,3) ,A,B,C,D,PRCT, COUNT 
INTEGER  I , J ,K ,L 


OPENING  FILES  FOR  READING  AND  WRITING 
CALL  EXCMS  (  '  FILEDEF  7  DISK  SUCVECT  DATA  Al'  ) 

CALL  EXCMS  (  ’  FILEDEF  8  DISK  FIRST  DATA  Al’  ) 

CALL  EXCMS  (  '  FILEDEF  9  DISK  SECOND  DATA  Al'  ) 

CALL  EXCMS  (  '  FILEDEF  15  DISK  BONUS  DATA  Al'  ) 


* 


10 


INITIALIZATION 
COUNT  =1.0 
BFI  =  0.0 
BSI  =  0.0 
BLINE  =  0.9250000 
DO  10  1=1,4 

SIGN(I)=.TRUE. 
FLAG( I)= .TRUE . 
CONTINUE 


*  READING  FIRST  INSPECTION  LOWER  CONFIDENCE  BOUNDS  FROM  DATA  FILE 
DO  50  1=1,2 
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DO  40  J=1 , 2 

DO  30  K=1 , 3 

DO  20  L=1 ,2 

READ(8 ,*)  FILCB( I , J ,K ,L) 
20  CONTINUE 

30  CONTINUE 

40  CONTINUE 

50  CONTINUE 


*  READING  SECOND  INSPECTION  LOWER  CONFIDENCE  BOUNDS  FROM  DATA  FILE 
DO  90  1=1,3 

DO  80  J=1 ,3 

DO  70  K=1 ,5 

DO  60  L=1 , 3 

READ (9 ,*)  SELCB( I , J ,K ,L) 

60  CONTINUE 

70  CONTINUE 

80  CONTINUE 

90  CONTINUE 

*  READING  SUCCESS  PROBABILITY  OF  FIRM  IN  EACH  TEST 

WRITEC*,*) 'WRITE  THE  PROBABILITY  OF  SUCCESS  IN  MANUFACTURER  TEST' 
READ  (*,*)  PSIM 

WRITEC*,*) 'WRITE  THE  PROBABILITY  OF  SUCCESS  IN  TEMPERATURE  AND 
+HUMIDITY  TEST' 

READ  (*,*)  PSITH 

WRITE(*,*) 'WRITE  THE  PROBABILITY  OF  SUCCESS  IN  VIBRATION  TEST' 
READ  (*,*)  PSIV 

WRITEC*,*) 'WRITE  THE  PROBABILITY  OF  SUCCESS  IN  ALTITUDE  TEST' 

READ  C*,*)  PSIA 

*  **1;+irkiT***4rto***irkit**irki;i;-kiT******\\-*1r.W".\-\-**-kiz*it-lc-k****1r*1z***-k****** 

*  READING  NUMBER  OF  SUCCESFUL  ITEMS 

*  FOR  FIRST  INSPECTION  AND  SECOND  INSPECTION 
DO  100  1=1,2000 
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READ(7 ,*)  SOM ( I ) , SOTH ( I ) , SOV ( I ) , SOA ( I ) 


100  CONTINUE 


JN.  JU ^LJLJUjL«jLJU^LA  JU ^ J* -L X X y^i «L XJ.. X J« «L *<*t"/t*^>1 ^  ^  "}\'}\  'A  M  ^  ^ 1 


INSPECTIONS  BEGIN 
DO  120  N=l, 2000,2 


*  SUBSCRIPT  DEFINITION  FOR  FIRST  INSPECTION 

A=INT( 20 . O-SOM(N) )+l 
B=INT(20 . O-SOTH(N) )+l 
C=INT( 32 . O-SOV(N) )+l 
D=INT(20 . O-SOA(N) )+l 

*  CHECK  FOR  ACCEPTANCE  OF  FIRST  OFFERED  LOT 

IF  (A.GT.2)  THEN 
SIGN(1)=. FALSE. 

END  IF 

IF  (B.GT.2)  THEN 
SIGN(2)=. FALSE. 

END  IF 

IF  (C.GT.3)  THEN 
SIGN(3)=.FALSE. 

END  IF 

IF  (D.GT.2)  THEN 
SIGN(4)=. FALSE. 

END  IF 

I F ( . NOT . ( S IGN ( 1 ) . AND . S IGN ( 2 ) . AND . S IGN ( 3 ) . AND . S IGN (  4  )  ) )  THEN 
SIGN( 1)=.TRUE. 

SIGN(2)=.TRUE. 

SIGN(3)=.TRUE. 

SIGN(4)= . TRUE . 

GO  TO  112 
END  IF 

*  'kicizJf'kJcicIrirk'iTii  Jrtcjctrk  A  /."A  A  A  A  /;  A  A  A  /»■ 

*  AFTER  FIRST  INSPECTION  DETERMINATION  OF  LCB  OF  FIRM 
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LCB=FILCB(A , B ,C ,D) 


IS  LOT  WORTH  WHILE  FOR  GETTING  BONUS  ? 
IF  (  LCB.GT.BLINE  )  THEN 
BFI=BFI+1 . 0 
GO  TO  110 
END  IF 


SUBSCRIPT  DETERMINATION  FOR  SECOND  INSPECTION 
A=A+INT(20.0-SOM(N+1)) 

B=B+INT( 20 . 0-SOTH(N+1) ) 

C=C+INT(32 . 0-SOV(N+1) ) 

D=D+INT(20 ,0-SOA(N+1) ) 

CHECKING  RESULTS  OF  SECOND  TEST  SERIES  ABOUT  ACCEPTANCE 
IF  (A.GT.3)  THEN 
FLAG(1)=. FALSE. 

END  IF 

IF  (B.GT.3)  THEN 
FLAG(2)=. FALSE. 

END  IF 

IF  (C.GT.5)  THEN 
FLAG(3)=.FALSE. 

END  IF 

IF  (D.GT.3)  THEN 
FLAG(4)=. FALSE. 

END  IF 

I F ( . NOT . ( FLAG ( 1 ) . AND . FLAG ( 2 ) . AND . FLAG ( 3 ) . AND . FLAG ( 4 ) ) )  THEN 
FLAG ( 1 ) = . TRUE . 

FLAG(2)=.TRUE . 

FLAG(3)=.TRUE. 

FLAG(4)=.TRUE . 

GO  TO  112 
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END  IF 


* 


* 


AFTER  SECOND  INSPECTION  DETERMINATION  OF  LCB 
LCB=SELCB(A ,B ,C ,D) 


*** 


CHECKING  FOR  BONUS  AFTER  SECOND  INSPECTION 
IF  (  LCB.GT.BLINE  )  THEN 
BSI=BSI+1 
GO  TO  112 
END  IF 


COUNTING  FOR  CHECKING  1000  REPLICATIONS 
110  COUNT  =  COUNT  +1.0 
A=0 . 0 


B=0 . 0 


C=0.0 


D=0 . 0 


* 


*  CHECKING  FOR  1000  REPLICATIONS 

IF  (COUNT. GT. 1000)  THEN 
GO  TO  130 
END  IF 

120  CONTINUE 

*  A-AAAA'A-AA'A-A-A'AAAA  A^nA-AAr-A-AArAAr-A  A  A"A  A  AAHlrArArAriHWIr*  A  A-A  A'AArAA'AA  A  A  A  AA.  A  A  A'A  AAArAAA 

*  PERCENTAGE  ESTIMATION  OF  GETING  BONUS  FOR  FIRM  A 
130  BTOT  =  BFI+BSI 

PRCT=BTOT/ 1000.0 

A  AAA**A***^HrtH:A-AA'A-AAA-A"AAAAAA  A  A  A  A  A-A  A  AAA-A  A  A  A'A  A  AAA-A-A-A-AA  A  AAAAAAA  A-AAAAA 

*  WRITING  RESULTS 
WRITE  (15,1) 

1  FORMAT  (// , 16X, '  BONUS  PLAN  SIMULATION  FOR  FIRM  A  ' ,2X) 

WRITE  (15,2) 

WRITE  (15,3) 
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F0RMAT(4X , '  FIRM  A  HAS  FOLLOWING  LONG  RUN  PROBABILITIES  IN  TESTS  ' 
+  ,2X) 

WRITE  (15,4) 

+  ./.2X) 

WRITE  (15,5)  PSIM,PSITH,PSIV,PSIA 

F0RMAT(4X, 'PROBABILITY  OF  SUCCESS  IN  MANUFACTURER  TEST  IS',2X, 

+F8. 6, //,4X, 'PROBABILITY  OF  SUCCESS  IN  TEMP.  AND  HUM.  TEST  IS',2X, 
+F8. 6, //,4X, 'PROBABILITY  OF  SUCCESS  IN  VIBRATION  IS’,2X, 

+F8. 6, //,4X, 'PROBABILITY  OF  SUCCESS  IN  ALTITUDE  IS ' ,F8 . 6 ,2X , //) 
WRITE  (15,6)  BFI 

F0RMAT(4X, 'FIRM  A  GAT  BONUS  AFTER  FIRST  INSPECTION  ' ,2X,F6. 1 ,2X, 

+ ' TIMES ' , 2X) 

WRITE  (15,7) 

+******', 2X,/) 

WRITE  (15,8)  BSI 

FORMAT(4X, 'FIRM  A  GAT  BONUS  AFTER  SECOND  INSPECTION  ' , 2X,F6 . 1 , 2X , 
+' TIMES* , 2X) 

WRITE  (15,9) 


9  FORMAT  (4X,  '***' 

+******' ,2X,/) 

WRITE  (15,11)  BTOT 

11  F0RMAT(4X, 'TOTALLY  FIRM  A  GAT  BONUS  IN  1000  REPLICATIONS  2X, 
+F5 . 1 ,2X, 'TIMES' ,2X) 

WRITE  (15,12) 

+******' ,2X,/) 

WRITE  (15,13)  PRCT 

13  F0RMAT(4X, 'GETTING  BONUS  PERCENTAGE  OF  FIRM  A  IS ' , IX ,F6 . 3 , 2X) 
WRITE  (15,14) 

14  F0RMAT(4X,'  ****-k*’t'"i>m**************'/rt''*************'‘ 

+^****' ,2 X,///) 

WRITE  (15,15)  BLINE 
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15  FORMAT (4X, '  BONUS  LINE  FOR  FIRMS  IS' , 1X,F6.4,2X) 
WRITE  (15,16) 

+******’ ,2 X,///) 

STOP 

END 
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APPENDIX  I.  95  %  LCB'S  FOR  DSBS  (  EQUAL  PROBABILITIES  ) 


Table  18.  95  %  LCB'S  FOR  DOUBLE  SAMPLING  BONUS  SYSTEM 


FAILURE 

VECTOR 

95  %  LCB 

FAILURE 

VECTOR 

95  %  LCB 

(0000) 

1.0000000 

(  0  0  0  1  ) 

0.9099184 

(0002) 

0.8706521 

(0010) 

0.9375849 

(001  1  ) 

0.8808635 

(0012) 

0.8387057 

(0020) 

0.8968240 

(  0  0  2  1  ) 

0.8528913 

(0022) 

0.8141473 

m 

0.8717731 

(  0  0  3  1  ) 

0.8281843 

(  0  0  3  2  ) 

0.7881665 

(0040) 

0.8422214 

(  0  04  1  ) 

0.8021229 

(0042) 

0.7645337 

(0100) 

0.9099184 

(0101) 

0.8548708 

(0102) 

0.8224728 

(0110) 

0.8808635 

(0111) 

0.8368122 

(0112) 

0.8010584 

(0120) 

0.8548636 

(0121) 

0.8131180 

(0122) 

0.7779465 

(0130) 

0.8294836 

(0131) 

0.7867751 

(0132) 

0.7541937 

(0140) 

0.8001103 

(0141) 

0.7624319 

(0142) 

0.7305543 

■nmm 

0.8657608 

(  0  2  0  1  ) 

0.8224728 

(  0  2  0  2  ) 

0.7905909 

(02  10) 

0.8387058 

(0211) 

0.7996263 

(0212) 

0.7653871 

(0220) 

0.8154084 

(  02  2  1  ) 

0.7777317 

(0222) 

0.7414687 

(  0  2  3  0  ) 

0.7881665 

(  0  2  3  1  ) 

0.7541937 

(  0  2  3  2  ) 

0.7210189 

(0240) 

0.7624319 

(  0  2  4  1  ) 

0.7326533 

(  0  2  4  2  ) 

0.6996729 

(  1  0  0  0  ) 

0.9945652 

(10  0  1) 

0.9099184 

(  1  0  0  2  ) 

0.8706521 
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Table  19.  95  %  LCB'S  FOR  DOUBLE  SAMPLING  BONUS  SYSTEM 


FAILURE 

VECTOR 

95  %  LCB 

FAILURE 

VECTOR 

95  %  LCB 

(1010) 

0.9375849 

(10  11) 

0.8808635 

(1012) 

0.8369564 

(  1  0  2  0  ) 

0.8968240 

(10  2  1) 

0.8519021 

(1022) 

0.8111921 

(  1  0  3  0  ) 

0.8717731 

(10  3  1) 

0.8258852 

(  1  0  3  2  ) 

0.7881665 

(  1  0  4  0  ) 

0.8422214 

(10  4  1) 

0.8001103 

(  1  0  4  2  ) 

0.7631623 

(110  0) 

0.9099184 

(110  1) 

0.8548708 

(1102) 

0.8224728 

(1110) 

0.8808635 

(Mil) 

0.8368122 

(1112) 

0.8002886 

(1  12  0) 

0.8528913 

(112  1) 

0.8121204 

(1122) 

0.7768002 

(113  0) 

0.8258852 

(113  1) 

0.7862302 

(1132) 

0.7532073 

<  1  14  0) 

0.8001103 

(114  1) 

0.7616678 

(1142) 

0.7295684 

(  1  2  0  0  ) 

0.8657608 

(1201) 

0.8224728 

(  1  2  0  2  ) 

0.7872553 

(1210) 

0.8387058 

(1211) 

0.7987179 

(1212) 

0.7626535 

(  1  2  2  0  ) 

0.8130434 

(1221) 

0.7758973 

(  1  2  2  2  ) 

0.7406143 

(  1  2  3  0  ) 

0.7860773 

(1231) 

0.7538874 

(  1  2  3  2  ) 

0.7200024 

(  1  2  4  0  ) 

0.7624319 

(1241) 

0.7309083 

(  1  2  4  2  ) 

0.6994067 

(2000) 

0.9891304 

(  2  0  0  1  ) 

0.9099184 

(2002) 

0.8657608 

(2010) 

0.9324048 

(2011) 

0.8808550 

(2012) 

0.8341966 
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Table  20.  95  %  LCB'S  FOR  DOUBLE  SAMPLING  BONUS  SYSTEM 


FAILURE 

VECTOR 

95  %  LCB 

FAILURE 

VECTOR 

95  %  LCB 

(  2  0  2  0  ) 

0.8968240 

(  2  0  2  1  ) 

0.8499787 

(  2  0  2  2  ) 

0.8110244 

(  2  0  3  0  ) 

0.8717731 

(  2  0  3  1  ) 

0.8247706 

(2032) 

0.7845937 

(  2  0  4  C  ) 

0.8422214 

(  2  041) 

0.8001103 

(2042) 

0.7624320 

(2100) 

0.9099184 

(2101) 

0.8548708 

(2102) 

0.8224728 

(2110) 

0.8808550 

(2111) 

0.8362838 

(2112) 

0.8003395 

(2120) 

0.8519021 

(2121) 

0.8121204 

(2122) 

0.7750615 

(2130) 

0.8258852 

(2131) 

0.7846695 

(2132) 

0.7517914 

(2140) 

0.8001103 

(2141) 

0.7605843 

(2142) 

0.7284731 

| 

1 

0.8608695 

(  2  2  0  1  ) 

0.8222825 

(2202) 

0.7860733 

(22  10) 

0.8341965 

(2211) 

0.7965909 

(22  12) 

0.7626535 

(2220) 

0.8110244 

(  2  2  2  1  ) 

0.7750614 

(2222) 

0.7401020 

(2230) 

0.7845957 

(2231) 

0.7519901 

(  2  2  3  2  ) 

0.7188347 

(2240) 

0.7611454 

(  2  2  4  1  ) 

0.7295684 

(2242) 

0.6975686 
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APPENDIX  J.  DOUBLE  SAMPLING  BONUS  SYSTEM  WITH  EQUAL 

PROBABILITIES 


Table  21.  DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.825 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.825 

PROB.S 

BONUS  % 

PROB.S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.008 

0.9375 

0.030 

0.9450 

0.059 

0.9500 

0.079 

0.9600 

0.172 

0.9700 

0.351 

0.9750 

0.487 

0.9800 

0.632 

0.9850 

0.799 

0.9900 

mm 

0.9950 

mm 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  *  0.B25 


Figure  5.  Double  Sampling  Bonus  System  With  LCBFB  =  0.825 
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Bonus  percentages  are  tabulated  and  plotted  below  with  LCBFB  0.850 


Table  22.  DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.850 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.850 

PROB/S 

BONUS  % 

PROB/S 

BONUS  % 

PROB/S 

BONUS  % 

0.9200 

0.007 

0.9375 

0.025 

0.9450 

0.050 

0.9500 

0.063 

0.9600 

0.144 

0.9700 

0.293 

0.9750 

0.423 

0.9800 

0.565 

0.9850 

0.738 

0.9900 

0.887 

mm 

0.970 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.B50 


Figure  6.  Double  Sampling  Bonus  System  With  LCBFB  =  0.850 
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Bonus  percentages  are  tabulated  and  plotted  below  with  LCBFB  0.875 


Table  23.  DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.875 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.875 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.007 

0.9375 

0.023 

0.9450 

0.045 

0.9500 

0.057 

0.9600 

0.122 

0.9700 

0.236 

0.9750 

mm 

0.9800 

mm 

0.9850 

mm 

0.9900 

0.799 

0.9950 

0.941 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.875 


Figure  7.  Double  Sampling  Bonus  System  With  LCBFB  =  0.875 
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Bonus  percentages  are  tabulated  and  plotted  below  with  LCBFB  0.900 

Table  24.  DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.900 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.900 

PROB.'S 

BONUS  % 

PROB.S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.002 

0.9375 

0.010 

0.9450 

0.020 

0.9500 

0.027 

0.9600 

0.052 

0.9700 

0.138 

0.9750 

0.212 

0.9800 

0.319 

0.9850 

0.486 

0.9900 

0.671 

0.9950 

0.872 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.900 


Figure  8.  Double  Sampling  Bonus  System  With  LCBFB  =  0.900 
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Bonus  percentages  are  tabulated  and  ploted  below  with  LCBFB  0.950 


Table  25.  DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.950 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.950 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.002 

0.9375 

0.010 

0.9450 

MEM 

0.9500 

0.022 

0.9600 

0.045 

mm 

0.107 

0.9750 

0.163 

0.9800 

0.237 

0.9850 

0.362 

0.9900 

0.502 

0.9950 

0.695 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.950 


Figure  9.  Double  Sampling  Bonus  System  With  LCBFB  =  0.950 


140 


Bonus  percentages  are  tabulated  and  ploted  below  with  LCBFB  0.999 


Table  26.  DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.999 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.999 

PROB/S 

BONUS  % 

PROB/S 

BONUS  % 

PROB/S 

BONUS  % 

0.9200 

0.002 

0.9375 

0.007 

0.9450 

0.010 

0.9500 

0.0 11 

0.9600 

0.028 

0.9700 

0.067 

0.9750 

0.106 

0.9800 

0.166 

0.9850 

0.274 

0.9900 

0.407 

0.9950 

0.626 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.999 


Figure  10.  Double  Sampling  Bonus  System  With  LCBFB  =  0.999 
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APPENDIX  K.  95  %  LCB'S  FOR  DSBS  (DIFFERENT  PROBABILITIES) 

Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.825 


Table  27.  DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.825 


0.950 

0.975 

0.990 

■SB 

MEM 

0.079 

0.364 

0.632 

0.689 

0.975 

0.099 

0.487 

mam 

0.882 

0.990 

0.109 

0.550 

0.920 

0.968 

mam 

0.110 

0.558 

0.934 

MEM 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  «  0.825 


Figure  11.  Double  Sampling  Bonus  System  With  LCBFB  =  0.825 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.850 


Table  28.  DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.850 


0.950 

0.975 

0.995 

0.950 

0.063 

0.320 

0.612 

0.680 

0.975 

mm 

0.423 

0.802 

0.871 

0.990 

0.083 

0.476 

0.887 

0.957 

0.995 

0.084 

0.483 

0.900 

0.970 

DOUBLE  SAMPLING  BONUS  SYSTEM 

LCBFB  =  0.850 


Figure  12.  Double  Sampling  Bonus  System  With  LCBFB  =  0.850 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.875 


Table  29.  DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.875 


0.950 

0.975 

0.990 

0.995 

0.950 

0.057 

0.277 

0.664 

0.975 

0.067 

0.353 

0.721 

0.846 

0.990 

0.072 

0.388 

0.799 

0.929 

0.995 

0.072 

0.392 

0.809 

0.941 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  =  0.875 


Figure  13.  Double  Sampling  Bonus  System  With  LCBFB  =  0.875 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.900 


Table  30.  DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.900 


0.950 

0.975 

0.990 

0.995 

0.950 

0.027 

0.164 

0.470 

0.605 

0.975 

0.031 

0.212 

0.603 

0.778 

0.990 

0.034 

0.235 

0.671 

0.860 

0.995 

0.034 

0.238 

0.679 

0.872 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  a  0.900 


Figure  14.  Double  Sampling  Bonus  System  With  LCBFB  =  0.900 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.950 


Table  31.  DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.950 


0.950 

0.970 

mm 

0.995 

0.950 

0.022 

mm 

0.369 

0.506 

0.975 

wem 

mm 

0.502 

0.686 

0.990 

0.029 

0.182 

0.502 

0.686 

0.995 

0.029 

0.185 

0.509 

0.695 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  -  0.950 


Figure  15.  Double  Sampling  Bonus  System  With  LCBFB  =  0.950 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.999 


Table  32.  DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.999 


0.950  0.970  0.990 


0.950  0.0  II  0.057 


0.975 


0.990 


0.995 


0.950 

0.970 

0.0 11 

0.057 

0.021 

0.106 

0.024 

0.151 

0.026 

0.161 

0.302 


0.407 


0.463 


0.995 


0.229 


0.406 


.570 


0.626 


DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  «  0.999 


. © 

'  . — 


.99 

prob.’S 


Figure  16.  Double  Sampling  Bonus  System  With  LCBFB  =  0.999 
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